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ABSTRACT 


The  C-17  test  and  evaluation  community  has  been  testing  different  aircraft 
formation  geometries  in  search  of  a  configuration  which  minimizes  paratrooper  encounter 
with  the  wake  vortices  of  upstream  aircraft.  This  thesis  develops  a  simulation  tool  that 
the  C-17  test  and  evaluation  community  can  utilize  as  an  advanced  risk  assessment  model 
to  use  on  proposed  formation  geometries  prior  to  live  testing.  The  model  is  developed 
under  the  architecture  of  object-oriented  simulation  using  MODSIM  III  and  parallels 
similar  efforts  by  the  Aerod3mamic  Decelerator  Technology  community  in  creating 
object-oriented  counterparts  to  already  developed  trajectory  models  of  various  degrees  of 
freedom. 

This  thesis  develops  the  paratrooper  object  portion  of  the  simulation  model  while 
the  Petry  thesis  (1997)  develops  the  C-17  aircraft  and  vortex  objects.  Once  integrated 
with  the  Petry  C-17  aircraft  and  vortex  objects,  and  after  verification  and  validation,  the 
simulation  model  is  applied  to  a  simplified  airborne  operation  scenario  using  the  mean 
distance  of  paratrooper  impact  location  to  assembly  areas  and  DZ  dispersal  distribution  as 
MOEs  for  different  aircraft  formation  geometries.  Lateral  separation  is  shown  to  have  the 
most  influence  on  both  MOEs,  while  trail  distance  has  minimal  effects.  For  the  airborne 
commander,  this  translates  into  operational  parameters  applicable  to  the  choice  of 
assembly  areas  and  formation  geometries.  Further  operational  parameters  of  any 
significance  are  gained  when  coupled  with  the  results  from  Petry  on  encounter  rates 
between  paratroopers  and  wake  vortices  where  trail  distance  has  a  significant  impact. 


XI 


C-17/PARATROOPER  RISK  ASSESSMENT  ANALYSIS 

1.  INTRODUCTION 

This  thesis’  purpose  is  to  provide  the  C-17  test  and  evaluation  community  with 
the  capability  to  assess  paratrooper  performance  during  C-17  drop  formations  in  both 
training  and  combat  environments.  Paratrooper  performance  takes  into  account  the  risk 
of  wake  vortex  encounter  and  ground  dispersal  patterns  on  the  drop  zone  (DZ).  Object- 
oriented  modeling  is  used  to  convert  current  static/deterministic  parachute/payload 
system  trajectory  models  of  any  degree  of  freedom  into  dynamic/stochastic  models 
through  the  development  of  a  class  of  parachute/payload  system  objects  which  are 
expandable  to  model  not  only  personnel  but  equipment  and  different  types  of  parachutes. 
The  immediate  impact  of  this  thesis  is  assessing  the  risk  of  C-17  formations  for  brigade- 
size  personnel  airborne  operations.  However,  the  parachute/payload  system  objects  can 
be  expanded  for  use  in  a  combat  modeling  environment. 

1.1  Background 

Paratrooper  interaction  with  wake  vortices  has  been  studied  by  the  US  Air  Force 
since  1987  conducting  the  first  tests  to  determine  the  effects  of  aircraft  wake  vortices  on 
parachute/payload  systems  using  static-line  deployed  parachutes  which  included 
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personnel  and  equipment  (Johnson  1988a:vii).  These  first  tests  investigate  the  effects  of 
C-130  and  C-141  aircraft  on  parachute/payload  systems,  while  follow-on  tests  investigate 
the  effects  from  C-5  wake  vortices.  Conclusions  from  both  studies  show  that  the  aircraft 
generate  significant  wake  vortices  in  size,  velocity  (strength) ,  and  duration.  Although  the 
wake  vortices  dissipate  with  time,  there  still  exists  the  danger  of  encounter  with 
parachute/payload  systems  in  formation  air  drops.  When  interacting  with  a 
parachute/payload  system,  these  wake  vortices  can  induce  potentially  hazardous 
conditions  to  include  parachute  collapse,  partial  deflation,  severe  oscillation,  increased 
rate  of  descent,  collision,  entanglement,  and  hard  landings  (Johnson  1988b;  13- 14). 
Recommendations  of  both  studies  indicate  the  need  to  develop  and  evaluate  formation 
tactics  for  large  scale  airborne  drop  activities,  along  with  safety  guidelines  to  avoid 
conditions  conducive  to  paratrooper/wake  vortex  encounter  and  interaction.  The  C-5 
study  also  found  a  need  to  begin  study  on  C-17  wake  vortex  characteristics  and 
parameters. 

With  the  advent  of  the  C-17  phasing  into  operational  roles  (as  a  C-141 
replacement)  it  has  encountered  several  operational  problems,  one  being  the  problem  of 
paratrooper/wake  vortex  encounter  when  used  as  a  jump  platform  for  US  Army  Airborne 
units.  During  operational  testing  of  airborne  activities  with  the  C-17  in  June  1995,  a 
paratrooper/wake  vortex  encounter  prompted  the  re-examination  of  formation  tactic 
geometry  (Blake  1996:1).  Due  to  the  size,  weight,  and  the  aerodynamics  of  the  C-17,  the 
wake  vortex  strength  is  greater  than  those  previously  encountered  with  the  C-130  or  the 
C-141  aircraft  (Natick  1996c:l).  As  with  the  previous  studies’  recommendations,  these 
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wake  vortices  are  key  considerations  in  developing  and  evaluating  formation  tactics  for 
mass  airborne  airdrop  activities  involving  paratroopers  and  their  equipment. 

In  developing  formation  tactics.  Air  Force  and  Army  planners  decided  that  the 
following  conditions  must  be  met:  (1)  total  wake  vortex  avoidance  must  be  assured 
where  the  wake  vortex  strength  presents  an  unacceptable  risk  to  the  paratrooper;  and,  (2) 
if  wake  vortex  avoidance  cannot  be  assured,  the  encounter  must  not  take  place  before  the 
wake  has  decayed  to  the  point  where  risks  associated  with  wake  vortex  strength  and 
probability  of  encounter  are  acceptable  (Natick  1996c:  1).  In  the  summer  of  1996,  a  series 
of  tests  were  conducted  at  the  drop  zones  located  at  Edwards  AFB,  CA  (EAFB)  in  an 
effort  to  find  the  probabilities  associated  with  wake  vortex  encounter  by  flying  in 
formations  most  conducive  to  paratrooper/wake  vortex  interaction.  A  total  of  57  test 
passes  were  made  with  a  total  of  672  dummies  dropped  from  both  C-141  and  C-17 
aircraft.  The  results  of  these  tests  were  carried  over  in  the  development  and  evaluation  of 
formation  tactics  at  Fort  Bragg,  NC  in  late  summer  of  1996.  Seeking  to  avoid 
paratrooper/wake  vortex  encounter,  these  tests  used  the  only  working  model  that  predicts 
paratrooper/wake  vortex  encounters  specific  to  the  C-17  (Blake  1996).  The  formation 
tactics  still  create  the  conditions  which  induced  paratrooper/wake  vortex  interaction, 
extending  the  testing  of  new  formation  tactics  until  the  conditions  set  forth  by  the  Air 
Force  and  Army  planner  can  be  met. 
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1.2  Problem  Statement 


Costs  of  test  and  evaluation  precludes  test  of  new  formation  tactics,  therefore  the 
need  for  a  modeling  tool  which  can  assess  the  probability  of  encounter  between 
paratroopers  and  wake  vortices  becomes  essential.  The  use  of  object-oriented  simulation 
in  modeling  airborne  airdrop  activities  can  predict  the  probability  of  paratrooper/wake 
vortex  encounter  under  different  formation  tactic  geometries.  Two  aspects  are  involved 
in  the  simulation  of  the  paratrooper/wake  vortex  interaction  system;  (1)  modeling  the 
wake  vortices  specific  to  the  C-17;  and,  (2)  modeling  a  parachute/payload  system’s 
trajectory  from  exit  to  impact. 

1,3  Organization  of  Research 

Research  of  work  already  accomplished  can  be  generalized  into  three  categories: 
(1)  research  and  development  of  parachute/payload  systems  for  use  in  space  related 
applications  (i.e.,  interplanetary  vehicle  atmospheric-entry  deceleration,  space  launch 
rocket  booster  recovery);  (2)  pure  complete  fluid  dynamics  of  parachute/payload  systems’ 
behavior  studies;  and,  (3)  research  and  development  of  parachute/payload  systems  for 
airdrop  applications.  Due  to  the  wide  range  of  behavioral  characteristics  modeled  in  the 
three  categories,  the  research  for  this  thesis  effort  concentrates  mainly  on  the  last  two 
categories. 

The  modeling  of  the  parachute/payload  system  and  its  surrounding  equations  of 
motion  when  influenced  by  gravitational  and  aerodynamic  forces  has  been  tackled  by 
engineers  for  years.  The  wide  array  of  approaches  differ  in  levels  of  approximation. 
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numbers  of  degrees  of  freedom,  and  the  numbers  of  independent  bodies  in  the  system 
(Maydew  1991;  121).  The  primary  research  focuses  on  the  areas  of  parachute/payload 
system  dynamics  and  parachute  trajectory  models,  and  searches  for  models  with  various 
degrees  of  freedom  to  achieve  a  specific  degree  of  accuracy.  Secondary  objectives 
include  finding  studies  on  the  interaction  of  parachute/payload  systems  and  aircraft  wake 
vortices,  and  mathematical  models  defining  the  dynamics  of  the  parachute/payload 
system  in  all  stages  of  descent. 

The  scope  of  this  thesis  effort  revolves  around  modeling  a  class  of 
parachute/payload  systems’  trajectory  from  aircraft  exit  to  impact  using  MODSIM  III,  an 
object-oriented  simulation  language.  This  thesis  effort  parallels  research  efforts  in  the 
aeronautical  community  as  evidenced  by  object-oriented  simulation  of  parachute 
trajectories  emerging  at  the  forefront  of  a  new  wave  of  modeling  techniques  at  the  14 
Aimual  AIAA  Aerodynamic  Decelerator  Systems  Technology  Conference  for  1997.  A 
defining  characteristic  of  this  thesis  effort  is  converting  well-defined  static/deterministic 
models  into  dynamic/stochastic  objects  to  provide  a  better  method  of  assessing  risk  of 
parachute/payload  system/wake  vortex  encounter,  and  assessing  groimd  dispersal  patterns 
on  the  drop  zone  using  different  C-17  airdrop  formation  geometries. 

The  main  premise  driving  this  thesis  effort  is  that  there  are  stochastic  aspects 
which  play  an  integral  part  of  the  parachute/payload  system  trajectory  from  the  time  the 
system  exits  the  aircraft  through  open-canopy  descent  propagating  all  the  way  to  impact 
dispersed  around  a  target  drop  point.  To  simplify  initial  development  of  the 
parachute/payload  system  object,  the  system  is  assumed  to  be  a  point-mass  system  with 
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two  degrees-of-freedom  (2-DOF).  After  successful  completion  of  the  2-DOF  point-mass 
object,  a  rigid,  single-body  system  using  higher  DOF  replaces  the  2-DOF  model  creating 
a  6-DOF  model.  A  rigid,  single-body  system  assumes  that  the  parachute  axis  of 
symmetry  remains  aligned  with  the  longitudinal  axis  of  the  payload  (Maydew  1991:122). 


/ 


/  .CJ 

J  ^ 


Figure  1.  Rigid,  Single-Body  System 

A  further  assumption  is  that  the  payload  is  unable  to  steer  the  parachute  system.  The  end 
product  will  be  a  paratrooper  object  with  6-DOF  modeled  by  an  unconscious  rigid,  single¬ 
body  system. 

1.4  Thesis  Objectives 

MODSIM  III,  an  object-oriented  simulation  language,  provides  the  architecture 
with  which  to  model  the  parachute/payload  system.  Chapter  II  provides  the  information 
necessary  to  model  a  parachute/payload  system,  and  introduces  trajectory  models  of 
various  DOF.  Chapter  III  provides  in  detail  the  modeling  of  the  parachute/payload 
system  trajectory  using  MODSIM  III. 
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The  simulation  itself  is  developed  in  three  phases.  The  first  phase  develops  a 
single  parachute/payload  system  object  in  MODSIM  III  from  an  already  validated  and 
academically  accepted  deterministic  model.  The  second  phase  modifies  the 
parachute/payload  system  object  into  a  paratrooper  object,  and  validates  the  results  with 
the  US  Army  Natick  Research,  Development,  and  Engineering  Center.  The  final  phase 
integrates  the  paratrooper  objects  with  the  Petry  (1997)  C-17  aircraft  and  wake-vortex 
objects,  using  both  discrete  and  continuous  simulation,  to  simulate  mass  airborne  drop 
formations;  and,  measures  dispersal  patterns,  with  no  interaction  from  wake  vortices, 
using  conventional  X/Y  coordinates. 

Chapter  IV  presents  the  results  of  the  simulation  when  used  with  a  simplified 
airborne  drop  scenario:  (1)  the  mean  distance  from  company  assembly  areas;  and,  (2)  the 
dispersion  distribution  of  the  paratrooper  object  across  the  width  and  down  the  length  of 
the  DZ.  Chapter  V  summarizes  the  thesis  effort  and  suggests  follow-on  efforts  in 
expanding  the  capability  of  the  object-oriented  simulation. 


1-7 


2.  LITERATURE  REVIEW 


The  literature  review  covers  two  principal  areas  of  research  interest— the  fluid 
dynamics  of  parachute-payload  systems’  behavior,  and  the  development  and  modeling  of 
simulations  of  parachute-payload  systems’  trajectory  during  descent.  However,  an 
understanding  of  the  standard  operating  procedures  for  airborne  operations  is 
fundamental  to  the  understanding  of  the  need  for  reliable  paratrooper  trajectory  modeling. 

2.1  Airborne  Operations 

The  US  Army  airborne  commander  and  staff  employs  a  form  of  “top-down 
analysis,”  referred  to  as  the  backward  planning  sequence,  in  planning  successful  airborne 
operations  (82"^*  ABN  DIV  ASOP  1985:2-1).  The  commander  follows  basic  joint 
airborne  operations  planning  principles  to  allow  for  maximal  effectiveness  in  the 
execution  of  the  airborne  operation.  The  airborne  commander  has  a  wide  view  of  the 
operations,  taking  into  account  the  limitations  of  the  airborne  units,  the  appropriateness  of 
the  mission  for  airborne  forces,  and  the  responsibilities  of  the  joint  task  force  (JTF)  or 
theater  commander. 

There  are  four  basic  plans  in  the  backward  planning  sequence,  all  focused  toward 
the  most  effective  way  of  acquiring  the  objectives  (82"‘*  ABN  DIV  ASOP  1985:2-6). 
Each  plan  supports  the  one  preceding  it,  as  the  planning  process  backs  out  from  the 
objectives.  The  four  basic  plans  are: 
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Ground  Tactical  Plan. 


•  Landing  Plan 

•  Air  Movement  Plan 

•  Marshaling  Plan 

In  the  Ground  Tactical  Plan,  the  conduct  of  operations  in  the  objective  area  is 
carefully  orchestrated  with  the  specifics  of  how  airborne  forces  are  to  maneuver  on  the 
ground.  The  Landing  Plan  directly  supports  the  Ground  Tactical  Plan.  The  main  focus  of 
the  Landing  Plan  concentrates  on  the  dispersion  of  the  airborne  units  in  the  DZ  in  such  a 
way  as  to  be  effectively  positioned  to  move  efficiently  into  their  area  of  responsibility 
(AOR)  as  designated  in  the  Ground  Tactical  Plan.  The  Air  Movement  Plan  encompasses 
the  period  of  time  from  aircraft  loading  to  arrival  of  airborne  units  in  the  objective  area. 
Specific  flight  routes  are  assigned  with  the  order  of  flight  and  formation  geometry  in  the 
area  around  the  DZ  chosen  that  will  best  support  the  Landing  Plan.  Incorporated  into  the 
Air  Movement  Plan  are  specific  arrival  times  over  the  DZ.  The  Marshaling  Plan 
completes  the  inverse  planning  sequence.  The  Marshaling  Plan  sees  airborne  units 
assemble  and  complete  their  final  preparations  for  combat,  and  coordinates  the  airborne 
units’  movement  to  departure  airfields.  Also,  in  the  Marshaling  Plan  airborne  units  are 
cross-loaded  onto  the  aircraft  such  that  the  Landing  Plan  is  best  supported  given  the  Air 
Movement  Plan. 
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2.2  A  irdrop  Studies 

In  developing  airdrop  models,  particular  interest  in  the  backward  planning 
sequence  focuses  around  the  Landing  Plan  and  the  Air  Movement  Plan.  The  Landing 
Plan  concerns  itself,  in  part,  with  the  dispersion  of  airborne  units  on  the  DZ,  and  the  Air 
Movement  Plan  concerns  itself  with  formation  geometry.  A  fundamentally  sound  method 
of  modeling  the  movement  of  paratroopers  along  their  trajectories  will  give  an  accurate 
estimate  of  paratrooper  dispersion  along  the  DZ.  Dispersion  is  a  function,  in  part,  of 
aircraft  formation  geometry  during  airdrop. 

Martin  (1978)  addresses  the  issue  of  missed  distance  errors  experienced  during 
Canadian  Forces  Paradrop  Exercises.  Martin  develops  a  simple  model  which  meets  the 
following  self-imposed  requirements  -  the  model  reasonably  describes  the  activity  in 
question;  the  parameters  affecting  airdrop  performance  are  correctly  varied  over  the 
effective  ranges;  the  range  of  possible  outcomes  is  not  so  wide  as  to  impair  the  usefulness 
of  simulation;  and,  the  model  is  verified,  if  possible,  using  empirical  data  (Martin 
1978:223).  Martin,  instead  of  using  a  parachute-payload  trajectory  model,  uses  a 
commonly  used  approximation  method  called  the  Computed  Air  Release  Point  (CARP) 
system  of  release  to  model  the  impact  point  of  the  parachute-payload  system.  The 
simulation  involved  varying  those  input  parameters  that  are  considered  the  most 
important:  winds,  position  at  release  point,  and  ballistics  or  hesitation  error  used  in  CARP 
calculations.  Martin,  however,  does  not  consider  formation  geometry  since  the 
simulation  focused  only  on  individual  parachute-payload  systems  and  not  on  mass 
formations. 
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Further  studies  done  on  airdrop  operations  were  done  by  Johnson  (1988)  which 
focused  more  on  vortex  encounter,  rather  than  on  DZ  dispersion,  with  the  C-130,  C-141 
and  C-5  aircraft.  Blake  (1996)  went  further  into  paratrooper/wake  vortex  encounters  in 
an  effort  to  predict  the  relative  locations  of  paratroopers  and  wake  vortices  of  the  C-17 
aircraft  during  formation  airdrops.  The  Blake  Model  can  be  modified  to  take  into  account 
a  specific  formation  geometry;  however,  it  does  not  attempt  to  predict  dispersion  on  the 
DZ. 


2.3  Dynamics  of  a  Parachute-Payload  System 

Basic  deceleration  system  equations  can  be  derived  starting  with  Newton’s 
Second  Law  and  building  upon  its  fundamental  premise.  In  order  to  understand  the 
equations  behind  parachute-payload  system  trajectory  models,  a  basic  understanding  of 
the  equations  that  govern  trajectory  dynamics  is  needed.  Seaman’s  (1975)  memorandum 
develops  the  fundamental  equations  needed  to  understand  the  dynamics  of  a  decelerating 
system  such  as  a  parachute-payload  system.  Seaman  breaks  down  the  components  that 
act  upon  a  basic  decelerating  system  to  mathematically  model  its  trajectory.  Knacke 
(1992)  confirms  Seaman’s  derivation  process  by  developing  the  basic  system  of  equations 
describing  trajectory  dynamics. 

In  describing  a  more  specific  application  of  a  parachute-payload  system,  more 
equations  are  added  to  the  basic  system  of  equations  to  model  characteristics  associated 
with  the  application.  When  parachute-payload  systems  are  released  from  aircraft,  new 
factors  (or  degrees  of  freedom)  help  describe  the  system’s  trajectory,  such  as  translation 
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and  rotation  (Jones  1987:538).  Jones  builds  on  both  the  rigid  single-body  and  rigid  two- 
body  models  to  describe  the  oscillation  effects  which  translation  and  rotation  causes.  For 
an  even  more  specific  application,  Heinrich,  Noreen,  and  Saari  (1973)  describe  the 
characteristics  of  a  parachute-payload  system  using  a  static-line  deployed  T-10  personnel 
parachute. 

The  modeling  of  a  parachute-payload  system  trajectory  can  be  expanded  with  the 
hierarchical  break  down  of  the  components  that  make  up  the  parachute-payload  system. 
The  trajectory  that  a  parachute-payload  system  follows  can  be  broken  down  into  discrete 
stages.  The  stages  are  payload  free  flight  (before  the  parachute  is  deployed);  parachute 
deployment;  parachute  inflation;  system  deceleration  and  turnover;  and,  steady-state 
descent.  Within  each  stage,  different  analytical  techniques  can  be  used;  however, 
Maydew  (1991)  shows  that  trajectory  dynamics  is  applicable  throughout  all  phases 
(Table  1).  Maydew  defines  trajectory  dynamics  as  the  analytical  solution  of  a  set  of  time- 
dependent  differential  equations  that  describe  the  trajectory  of  a  system.  Maydew  lists 
several  references  that  discuss  the  equations  of  motion,  the  selection  of  an  axis  system, 
and  simplifying  assumptions.  He  furthers  lists  several  other  references  that  have  made 


Table  1 

Modeling  the  Stages  of  Parachute/Payload  System  Flight  (Maydew  1991:121) 


STAGES 

Modeling 

Technique 

Payload 
Free  Flight 

Decelerator 

Deployment 

Inflation 

System  Deceleration 
and  Turnover 

System 

Descent 

Trajectory  Dynamics 

• 

• 

• 

• 

• 

Elastic  Structure 

• 

• 

Dynamics 

Decelerator 

• 

• 

• 

• 

Steady-State 

Aerodynamics 

Decelerator  Unsteady 

• 

• 

Aerodynamics 
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significant  contributions  in  the  development  of  parachute-payload  system  trajectory 
models. 

Maydew  also  presents  several  models  that  describe  the  system  state  of  a 
parachute-payload  system’s  trajectory.  The  most  basic  trajectory  model  is  that  of  a  point 
mass  system  in  which  all  forces  act  upon  a  central  point  in  the  system  in  the  vertical 
plane.  This  is  considered  a  2-DOF  system  (Figure  2)  because  it  contains  only  horizontal 
and  vertical  components  to  describe  its  motion.  The  next  level  in  trajectory  modeling  is  a 
3-DOF  system  (Figure  3)  which  adds  an  oscillation  component  to  the  2-DOF  system  of 
equations,  although  it  is  still  limited  to  the  vertical  plane.  Within  the  family  of  3-DOF 
models,  Maydew  shows  that  the  system  itself  can  be  modeled  several  ways:  massless, 
rigid  single-body,  rigid  two-body,  and  elastic.  A  massless  decelerator  model  maintains 
the  alignment  of  the  decelerator  along  the  payload  velocity  vector,  and  the  only  force  that 
the  decelerator  produces  is  drag.  A  rigid  single-body  system  is  one  where  the  parachute 
axis  of  symmetry  remains  aligned  with  the  payload’s  longitudinal  axis;  and,  the  forces 
produced  by  the  decelerator  are  drag,  lift,  and  a  moment  in  which  all  are  assumed  to  act  in 
the  payload’s  center  of  gravity.  The  rigid  two-body  system  is  a  point  mass  pinned  to  a 


Figure  2.  2-DOF  System  Figure  3.  3-DOF  System 
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rigid  body  by  a  massless  rigid,  connection.  Finally,  the  elastic  model  has  mass  node 
equations  of  motion,  defining  the  decelerator  and  elastic  suspension  lines,  which  are 
solved  simultaneously  with  the  rigid  body  equations.  Doherr  (1992)  summarizes  the 
different  parachute  models  and  their  capabilities  (Table  2).  For  the  6-DOF  model,  y/,  6, 
and  (j)  describe  rotational  motion  of  the  system  (Figure  4).  The  9-DOF  model  adds  to  this 
the  variables  y/p,  dp,  and  (ftp,  which  describe  the  relational  and  rotational  motion  of  the 
payload  with  respect  to  the  canopy. 


<!> 


Figure  4.  6-DOF  System 
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Table  2 

Comparisons  of  Parachute  Models  (Doherr,  1992:6-4) 


Trajectory  Analysis 

Point  Mass 

Planar 

6-DOF 

9-DOF 

or  Ballistic 

Rigid  Body 

Rigid  Body 

Rigid  Body 

Degrees  of  Freedom  (DOF) 

2 

3 

6 

9 

Major  Variables 

X,  z 

X,  z,  0 

X,  y,  z, 

X,  y,  z, 

V,  6,  <t> 

V|/,  0,  (j), 

_ Vpi  6p»  *t>p _ 

PAYLOAD  INPUTS 

Mass 

• 

• 

• 

• 

Inertias 

• 

• 

• 

Cx 

• 

• 

• 

• 

Cy 

• 

• 

Cz 

• 

• 

• 

Cl 

• 

• 

Cm 

• 

• 

• 

Cn 

• 

• 

DECELERATOR  INPUTS 

Mass 

• 

Inertias 

• 

CdS  (drag  area) 

• 

• 

• 

• 

Cn  (normal) 

• 

C|  (roll) 

• 

• 

Xcp  (center  of  pressure) 

• 

ttii  (apparent  mass) 

• 

Coupling  Conditions 


OUTPUT 


Deceleration 

• 

• 

• 

• 

Velocity 

• 

• 

• 

• 

Down-Range 

• 

• 

• 

• 

Off-Range 

• 

• 

Altitude 

• 

• 

• 

• 

Heading  Angle 

• 

• 

Pitch  Angle 

• 

• 

• 

Roll  Angle 

• 

• 

Relative  Angles 

• 

In  communicating  trajectory-modeling  methods  in  the  airborne  test  and  evaluation 
community,  the  stages  of  system  flight  are  commonly  referred  to  as  zones  (Figure  1). 
Zone  1  is  defined  from  the  time  of  exit  to  canopy  inflation  just  after  first  vertical.  This  is 
normally  130  below  exit  for  the  T- IOC.  Zone  2  is  the  region  from  130  feet  below  exit  to 
450  feet  above  ground  level  (AGL).  In  combat  jumps,  Zone  2  is  non-existent.  Zone  3  is 
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the  region  from  450  AGL  to  ground.  In  combat  jumps.  Zone  3  is  defined  as  the  region 
from  130  feet  below  exit  to  ground  level  (Natick  1996). 


EXIT 


Figure  5.  Parachute  Behavior  Zone  Regions  (Natick  1996) 

2.4  Current  Models 

The  Blake  Model  (1996)  uses  a  simple  point-mass  system  with  3-DOF,  one  each 
for  the  vertical,  horizontal,  and  lateral  components  of  motion.  Although  the  model  uses 
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3 -DOF,  both  the  horizontal  and  lateral  components  are  simplified  by  using  the  same 
components  of  cross  wind  as  the  defining  relation. 

Another  3 -DOF  model  (Bermey  1996)  uses  the  rigid  two-body  model  that 
oscillates  in  a  two-dimensional  wind  profile.  This  model  uses  the  basic  system  of 
equations  developed  from  Newton’s  Laws  of  Motion  to  estimate  the  response  of  a 
parachute-payload  system  when  encountering  a  changing  wind  profile  such  as  a  wake 
vortex  approximation. 

Natick  RDEC  has  developed  another  3-DOF  model  for  a  different  purpose 
(Wallace  1996).  The  motivation  in  this  model  is  to  show  an  inherent  Dispersion  Error 
Probability  (DEP)  for  a  non-gliding  (unconscious)  parachute-payload  system.  The  model 
uses  a  point-mass  system  with  random  oscillation  affecting  both  lateral  and  horizontal 
components  of  the  trajectory.  Factors  that  affect  the  DEP  are  identified  and  are  taken  into 
account  if  considered  to  possess  a  significant  affect  on  DEP. 

A  6-DOF  model  increases  the  resolution  of  the  trajectory  and  can  account  for  the 
position  components  of  the  trajectory  which  are  the  vertical,  horizontal,  and  lateral 
movements  (with  corresponding  axes)  of  the  parachute-payload  system,  as  well  as  three 
added  DOF  for  rotation  about  each  axes.  One  such  model,  developed  by  Tory  and  Ayres 
(1977),  uses  a  rigid  two-body  system  and  models  the  trajectory  of  the  6-DOF  parachute- 
payload  system  after  full  canopy  inflation. 

Models  that  use  high  DOF  are  more  useful  when  analysis  of  the  relative  motion 
between  the  parachute  and  the  payload  needs  to  be  considered  (Doherr  1992:774). 
Doherr  (1992)  develops  a  more  complex  model  having  9-DOF,  which  uses  a  non-rigid 
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two-body  system  where  two  masses,  the  parachute  and  the  payload,  are  cormected  by  a 
joint.  There  are  6-DOF  associated  with  the  payload  (position  and  attitude)  and  3 -DOF 
associated  with  the  parachute  (attitude  only). 

This  is  by  no  means  an  exhaustive  list  of  models  available,  but  is  one  that  are 
useful  for  the  level  of  resolution  this  thesis  effort  aspires  toward.  In  a  complete  fluid 
dynamic  (CFD)  study  of  parachute  structural  dynamics  in  a  close  surface  investigation, 
thousands  of  DOF  are  required  due  to  the  ever-increasing  complexities  involved,  and 
therefore  is  beyond  the  scope  of  this  effort. 


2-11 


3.  METHODOLOGY 


There  exists  a  need  to  develop  a  usable  airborne  risk  assessment  model  helpful  in 
both  the  Landing  Plan  and  the  Air  Movement  Plan  in  the  backwards  plarming  sequence  of 
airborne  activities,  and  addressing  both  paratrooper/wake  vortex  encounter  and  DZ 
dispersion.  Given  the  current  trial  and  error  method  of  formation  testing  in  the  C-17 
airdrop  test  and  evaluation  community,  aided  by  the  Blake  (Wright  Lab)  Model,  the 
problem  is  to  develop  a  stochastic  paratrooper  object  which  predicts  paratrooper 
trajectories  to  be  used  in  conjunction  with  the  Petry  C-17  aircraft  and  vortex  objects.  The 
technique  uses  object-oriented  simulation  subject  to  similar  self-imposed,  but  necessary, 
requirements  of  Martin’s  (1978)  model. 

3.1  Solution  Technique 

Exhaustive  testing  of  possible  formation  geometries  only  provides  data  for  those 
specific  cases.  An  alternative  to  the  trial  and  error  approach  that  the  airborne  test  and 
evaluation  community  uses  is  to  approximate  the  complex  system  environment  of 
simultaneous  activities  and  events  involved  in  airdrop  operations  using  simulation 
methods.  Smart  and  effective  testing  methods  exist  with  the  use  of  experimental  design 
analysis  (EDA)  and  response  surface  methodology  (RSM);  but,  at  the  same  time 
questions  arise  as  to  whether  or  not  the  formation  geometry  being  tested  is  best  suited  for 
the  conditions  which  exist.  Given  the  cost  of  testing,  repetition  of  sample  runs  may  not 
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be  feasible.  Using  a  simulation  model,  repetitive  testing  can  be  performed  under  a  varied 
range  of  controllable  parameters  to  determine  feasible  configurations  prior  to  live  testing. 

Different  forms  of  simulation  methods  exist,  ranging  from  spreadsheets  to  highly 
specialized  languages  specifically  catered  for  simulation.  MODSIM  III  is  one  such 
specialized  simulation  language  in  that  it  is  object-oriented.  Objects  are  self-contained 
data  structures  that  have  their  own  methods  (Banks  et  al.  1996:  128).  They  can  be  looked 
upon  as  the  building  blocks  of  a  system,  behaving  independently  of  other  building  blocks 
yet  interacting  with  them.  When  taken  as  a  whole,  object-oriented  programming  provides 
a  more  natural  way  to  approach  the  building  of  a  system  simulation  by  dynamically 
adding  more  instances  of  an  object  as  they  are  needed.  Object-oriented  programming 
facilitates  multiple,  concurrent  instances  of  similar  objects. 

In  airborne  operations,  several  entities  exist  concurrently  with  each  entity 
exhibiting  specific  behaviors.  There  are  multiple  aircraft  in  formation,  dropping  multiple 
paratroopers  and  generating  wake  vortices,  all  in  the  same  airspace.  Although  the 
aircraft,  vortices,  and  paratroopers  behave  independently  of  each  other,  they  continually 
interact  with  one  another— the  aircraft  fly  in  a  formation,  the  paratroopers  exit  the  aircraft, 
and  the  wake  vortices  may  or  may  not  disturb  the  paratroopers’  trajectories.  In  short,  the 
type  of  system  environment  involved  with  airborne  operations  lends  itself  well  to  object- 
oriented  simulation. 
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3.2  Implementation 

US  Army  Natick  Research  Development  and  Engineering  Center  remains  at  the 
forefront  of  aerodynamic  decelerator  systems  modeling  specific  to  US  Army  needs.  From 
the  outset,  Natick  has  been  involved  in  the  development  of  the  paratrooper  object.  With 
Natick,  the  decision  to  pursue  a  6-DOF  model  consistent  with  Table  2  is  substantiated 
because  of  the  need  for  the  information  that  a  6-DOF  model  provides  over  a  3 -DOF 
model  and  because  a  9-DOF  model  would  provide  no  added  benefit  (Benney  1996).  The 
Purvis  6-DOF  Model  (1987),  written  in  FORTRAN,  is  a  simple  yet  effective  model 
providing  a  straight  forward  approach  to  trajectory  propagation  by  using  a  second-order 
Euler  integration  scheme  (Doherr  1992:6-1).  The  development  of  a  paratrooper  object 
takes  place  in  three  stages.  First,  the  original  Purvis  FORTRAN  Model  code  is  translated 
into  an  object  using  the  object-oriented  simulation  language  MODSIM  III.  Once  a 
successful  translation  is  achieved,  the  object  is  modified  to  reflect  the  aerodynamic 
properties  of  a  combat  equipped  paratrooper  using  a  T-IOC  parachute.  The  last  stage 
integrates  the  paratrooper  objects  with  the  Petry  C-17  aircraft  and  wake  vortex  object  and 
incorporates  the  addition  of  stochastic  elements  in  the  trajectory  propagation  scheme. 

3.2.1  Translating  Into  MODSIM  HI 

The  Purvis  Model  is  organized  using  a  single  main  program,  four  subroutines,  and 
a  function  routine.  However,  it  can  be  looked  upon  as  performing  three  distinct  roles: 
input  of  parameters  specific  to  the  parachute-payload  system,  initialization  of  starting 
conditions,  and  trajectory  propagation.  Similar  roles  are  needed  in  order  to  perform  the 
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same  calculations  in  MODSIM  III.  An  initial  structure  of  four  MODSIM  III  modules  has 
been  used  to  properly  convert  the  Purvis  Model.  The  modules  are: 

•  Main  Module 

•  Global  Module 

•  Calculation  Module 

•  Jumper  Module 

The  Main  Module  is  a  simple  module  that  most  accurately  describes  the  activities 
taking  place  in  the  simulation.  In  the  Main  Module,  a  new  instance  of  a  jumper  object, 
emulating  the  same  type  of  parachute-payload  system  used  in  the  Purvis  Model,  is  created 
and  told  to  jump  at  the  start  of  the  simulation.  Once  the  jumper  object  impacts  ground 
level,  it  is  disposed. 

The  Global  Module  defines  and  initializes  those  constants  and  looping  variables 
that  are  used  throughout  the  simulation.  It  includes  the  method  initializeData  that  simply 
initializes  the  constants  used  in  trajectory  propagation.  The  Main  Module  imports  this 
method  and  implements  it  prior  to  creating  a  new  instance  of  a  jumper  object. 

The  Calculation  Module  defines  two  functional  procedures  that  are  used 
repetitively  during  trajectory  propagation.  The  first  procedure  calculates  the  gravitational 
effect  at  a  given  altitude.  The  other  procedure  calculates  the  density  and  speed  of  sound 
for  a  given  altitude  above  sea  level  (ASL). 

The  Jumper  Module  contains  the  bulk  of  the  MODSIM  III  code.  There  are  two 
methods  written  into  the  Jumper  Module:  Objinit  and  jump.  When  the  Main  Module 
creates  a  new  instance  of  a  jumper  object,  the  Objinit  Method  is  automatically  called 
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upon  to  initialize  all  the  parameters  specific  to  the  jumper  in  its  trajectory  propagation. 
The  jump  Method  houses  the  trajectory  propagation  loop.  Trajectory  information 
continues  to  be  calculated  until  ground  impact.  During  this  propagation,  the  jumper 
maintains  all  output  information  consistent  with  Table  2  and  displays  the  information  at 
specified  time  intervals. 

The  difficulty  in  translating  from  the  Purvis  FORTRAN  Model  to  a  MODSIM  III 
object  is  two-fold.  Since  a  basic  understanding  of  the  logic  flow  of  the  Purvis  Model 
must  be  discerned  from  the  FORTRAN  coding,  the  first  level  of  difficulty  is 
understanding  of  the  structural  flow  of  logic,  i.e.,  what  is  the  Purvis  Model  doing?  The 
second  level  of  difficulty  is  translating  the  FORTRAN  structure  into  a  MODSIM  III 
architecture.  Both  difficulties  have  been  overcome  with  an  end  result  of  a  MODSIM  III 
object  that  correctly  emulates  the  FORTRAN  parachute-payload  system  of  the  Purvis 
Model. 

3.2.2  Modification  Into  Paratrooper  Objects 

Modifying  the  parachute-payload  system  of  the  Purvis  Model  is  a  matter  of 
modifying  the  input  parameters  to  correctly  represent  the  aerodynamic  properties  of  a 
combat  equipped  paratrooper  using  a  T-IOC.  Under  the  guidance  of  Natick,  the 
following  input  parameters  have  been  changed  in  order  to  model  a  paratrooper’s 
trajectory:  weight,  payload  (forebody)  center  of  gravity  (eg),  payload  length,  inertia  in  the 
x-coordinate  (roll),  inertia  in  the  y-coordinate  (pitch),  inertia  in  the  z-coordinate  (yaw), 
and  effective  T-IOC  inflation  area.  The  Purvis  Model  has  an  example  that  changes  all 
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these  aerodynamic  properties  in  its  input  routine  to  show  that  the  model  can  be  used  with 
different  parachute-payload  systems.  The  modifications  needed  in  the  development  of  a 
paratrooper  object  follow  closely  to  Purvis’  own  changes. 

Several  assumptions  are  made  regarding  the  physical  properties  of  the 
paratrooper-parachute  system  in  order  to  find  its  aerodynamic  properties.  The 
paratrooper  object  assumes  the  physical  geometry  of  a  cylinder,  with  the  T-IOC  assuming 
the  form  of  a  half-sphere.  This  is  a  typical  assumed  configuration  used  in  parachute- 
payload  systems  modeling.  There  is  also  the  conical  section  of  the  system  formed  by  the 
risers  connecting  the  paratrooper  to  the  parachute.  And  lastly,  as  the  parachute  inflates, 
there  is  added  mass  to  the  system  due  to  the  air  trapped  under  the  parachute  canopy.  The 
cylinder  is  assumed  to  be  six  feet  high  and  one  foot  wide.  Although  joined  at  a  point 
where  the  harness  meets  the  paratrooper,  the  system  is  considered  a  rigid  body.  In  other 
words,  the  joint  does  not  pivot,  and  the  payload  has  no  relative  motion  in  relation  to  the 
parachute.  Initially,  the  weight  of  a  360-pound  paratrooper  is  assumed  since  there  are 
known  parameters  for  such  a  paratrooper,  such  as  mean  descent  velocity  and  mean 
descent  time.  Moments  of  inertia  around  each  of  the  axes  are  found  easily  using  a  set  of 
equations  foimd  in  basic  physics  texts  (Benney  1996). 

In  finding  the  moments  of  inertia,  the  following  parameters  are  first  defined: 
a  -  radius  of  payload  cross-section 

dc  =  distance  from  canopy  eg  to  system  eg 
dp  =  distance  from  payload  eg  to  system  eg 
ds  =  distance  from  suspension  line  eg  to  system  eg 
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Ip  =  length  of  payload 
Ma  -  added  mass 
Mp  =  mass  of  payload 
r  =  radius  of  canopy 
p  =  air  density 
Wc  =  weight  of  canopy 
Wp  =  weight  of  payload 
Ws  =  weight  of  suspension  line 

First,  the  added  mass  of  an  inflated  canopy  is  calculated  using 

4  , 

M^=p-7vr  . 

The  payload  moment  of  inertia  about  the  x  and  y-axes  are  calculated  using 

and  about  the  z-axis  using 

I  =^Mpa\ 

Pzz  2 


Finally,  using  Ma  and  Ip,  the  system  moment  of  inertia  about  the  x  and  y-axes  are 


2  W  W  W 


and  the  z-axis  is 


2  W  2 
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The  last  change  to  the  p^ameters  involves  the  effective  inflation  area  of  the  T- 
lOC.  Lying  down  flat  the  coverage  area  of  the  T-IOC  is  953  (Natick  1996:1). 
However,  once  inflated,  the  effective  inflation  area,  assuming  the  inflated  T-IOC  becomes 
a  half-sphere,  is  approximately  690  ft^. 

3.2.3  Integration  and  Enhancements 

The  building  block  ease  of  MODSIM  III  allows  for  easy  integration  between 
different  objects.  The  paratrooper  object  combines  with  the  Petry  C-17  aircraft  and 
vortex  objects  to  form  a  cohesive  system  to  model  airborne  operations.  During  the 
integration  process,  the  need  for  further  modifications  becomes  apparent.  The  following 
new  modifications  are  made: 

•  management  of  positional  information  using  three  coordinate  systems 

•  separate  right  and  left  paratrooper  objects 

•  a  system  to  pass  relevant  information  between  the  C-17  aircraft  objects  and 
the  individual  paratrooper  objects 

•  a  greenLight  Method  which  initiates  the  paratrooper  object  jump  formation 

•  a  method  for  each  paratrooper  object  to  calculate  its  own  distance  from 
vortices  originating  from  upstream  aircraft 

•  the  introduction  of  stochastic  variables  into  the  paratrooper  objects 
Coordinate  Systems.  The  use  of  three  related  coordinate  systems  becomes  a 

convenient  way  of  managing  relative  positional  information  between  the  C-17  aircraft, 
vortices,  and  paratrooper  objects.  The  first  coordinate  system  has  its  origin  centered  on 
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the  lead  aircraft  of  the  lead  element,  is  referred  to  as  the  Aircraft  Coordinate  System 
(ACS),  and  is  continuously  moving.  Positive  direction  is  measured  aft  of  the  aircraft  for 
the  x-position,  starboard  for  the  y-position,  and  above  for  the  z-position.  The  Ground 
Coordinate  System  (GCS),  the  stationary  system,  has  its  origin  at  leading  edge  of  the  drop 
zone.  Positive  direction  is  measured  in  the  direction  of  flight  path  for  the  x-position,  to 
the  right  for  the  y-position,  and  above  for  the  z-position.  A  third  coordinate  system,  the 
Inertial  Coordinate  System  (ICS),  is  unique  to  the  paratrooper  objects.  Its  origin  is  at  the 
point  of  exit  for  the  paratrooper  object  and  remains  stationary  from  exit  to  impact. 
Positive  direction  is  measured  along  the  direction  of  flight  for  the  x-position,  starboard 
for  the  y-position,  and  downward  for  the  z-direction. 

Left  and  Right  Paratrooper  Objects.  The  addition  of  distinct  right  and  left 
paratrooper  objects  adds  to  the  resolution  of  the  simulation  by  modeling  the  exit  of 
paratroopers  from  either  the  right  or  left  rear  doors  of  the  C-17.  The  only  difference 
between  a  right  and  a  left  paratrooper  object  is  its  point  of  origin. 

Communication  of  Information.  Communication  between  the  aircraft  and  the 
paratrooper  objects  is  essential  for  several  activities  that  occur  in  the  simulation.  The  first 
instance  of  shared  information  is  in  the  greenLight  Method  of  the  C-17  aircraft  objects. 
The  greenLight  Method  starts  the  stick  of  paratrooper  objects  exiting  the  aircraft.  In  this 
method,  the  C-17  aircraft  object  passes  on  its  positional  information  to  the  paratrooper 
objects.  Once  free  of  the  aircraft,  the  paratrooper  objects  act  independently  of  any  other 
object;  however,  several  pieces  of  information  need  to  be  continuously  passed  into  the 
paratrooper  objects,  since  they  change  concurrently  with  the  paratrooper  objects’ 
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trajectories.  MODSIM  III  allows  for  easy  exchange  of  information  between  different 
objects  as  long  as  the  requesting  object  specifies  which  object  to  import  the  information 
from.  Table  3  lists  the  categories  of  information  which  are  continuously  changing  and  are 
continually  used  by  each  paratrooper  object. 

Table  3 

Location  of  Imported  Information 

"Imported  from 

Vortex  Control  Module 
Vortex  Control  Module 
Global  Module 
Vortex  Module 

Green  Light.  In  airborne  operations,  aircraft  are  flying  in  a  specified  formation. 
Green  Light  is  called  when  the  CARP  is  encountered.  The  CARP  in  the  simulation  is 
similar  to  the  CARP  in  the  model  done  by  Martin  in  that  it  is  that  point  in  airspace  where 
paratrooper  objects  begin  their  exit  from  the  aircraft.  However,  it  differs  in  that  it  is  not 
calculated  from  a  planned  impact  point  on  the  DZ,  but  rather,  it  is  used  as  a  starting  point 
for  the  simulation.  Once  Green  Light  is  reached,  paratrooper  objects  begin  to  exit  the 
aircraft  until  the  last  paratrooper  object  in  the  stick  has  departed.  Paratrooper  objects  exit 
in  a  static  inter-departure  time  of  0.5  seconds. 

Vortex  Polling.  Management  of  positional  information  allows  for  each 
paratrooper  object  to  calculate  its  position  in  all  three  coordinate  systems  in  order  to  keep 
track  of  relative  positions.  Relative  positioning  is  used  in  the  determination  of 
paratrooper  object  distances  from  known  vortex  positions  called  missed  distance.  This  is 
accomplished  in  the  pollVortices  Method.  The  linear  algebra  method  of  vector  projection 
is  used  to  find  the  orthogonal  distance  between  the  paratrooper  object  and  the  vortex 


INFORMATION 
Exit  Positions 

Position  of  Lead  Aircraft  in  Lead  Element 
Variable  Winds 
Vortex  Positions 


object’s  core.  The  core  positions  of  vortex  objects  are  defined  as  points  in  space  at  100- 
foot  intervals,  originating  100  feet  behind  its  generating  aircraft.  Since  the  vortex  objects 
are  defined  under  the  ACS,  each  paratrooper  object  must  translate  its  position  from  the 
ICS  to  the  ACS  to  begin  a  search  along  the  vortex  body  for  its  closest  orthogonal 
distance.  Before  calculating  orthogonal  distances,  the  paratrooper  object  searches  for  that 
position  along  the  vortex  body  where  it  is  within  a  tolerance  box;  i.e.,  a  region  that  is 
reasonably  close  to  the  vortex.  Once  within  this  tolerance  box,  the  paratrooper  object 
calculates  the  orthogonal  distance  by  first  defining  two  vectors  having  the  same  origin. 
The  first  vector  is  defined  from  a  vortex  object  core  position  to  the  paratrooper  object 
position  V. .  The  second  vector  is  defined  from  the  same  vortex  core  position  to  the 

subsequent  vortex  core  position  v, .  In  Figure  6,  projection  of  Vj  onto  v,  is  calculated 
using  the  following  equation  (Strang  1986:  147) 


where  is  the  projection,  and  the  orthogonal  distance  is  Ivy-vJ.  If  this  orthogonal 

distance  is  within  the  effective  vortex  radius  (defined  by  a  critical  vortex  strength,  or 
swirl  velocity)  then  an  encounter  is  recorded.  Major  and  minor  encounters  are  not 
distinguished.  The  paratrooper  object  makes  these  calculations  for  every  upstream  vortex 
generated. 
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Figure  6.  Orthogonal  Distance  of  Paratrooper  from  Vortex  Core 


Inclusion  of  Stochastic  Elements.  Without  any  form  of  stochastic  elements 
present  in  the  model,  the  vortex  encoimter  rate  and  the  DZ  dispersion  can  be  calculated 
deterministically.  The  trajectories  of  each  paratrooper  object  mirrors  exactly  the 
trajectory  of  any  other  paratrooper  object  resulting  in  either  all  or  no  paratrooper/wake 
vortex  encounters.  The  addition  of  stochastic  behavior  enhances  the  simulation’s 
usefulness  by  modeling  their  random  behavior.  The  variables  to  randomize  are  chosen 
such  that  the  model  reflects  actual  occurrences  in  an  airborne  operation.  Although  any 
numbers  of  elements  are  candidates  for  stochastic  modeling,  this  thesis  limits  the  choices 
to  paratrooper  weights  and  T-IOC  glide.  In  keeping  with  resolution  and  capabilities  of 
the  Purvis  Model,  these  two  elements  are  easily  modeled  without  changing  the 
fundamental  aerodynamic  methodology  used.  (With  the  integration  of  aircraft  and  vortex 
objects,  winds  are  a  stochastic  element  that  affects  trajectory  propagation.)  From  the 
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weights  of  paratroopers  used  in  D-bag  clearance  testing  in  March  1996,  a  normal 
distribution  is  fitted  to  model  paratroopers  with  a  mean  weight  of  247  pounds  and  a 
standard  deviation  of  24.35  pounds. 

A  harder  problem  is  the  random  glide  inherent  in  the  T-IOC  when  no  wind  is 
present  and  no  oscillation  is  being  experienced.  Natick  recognizes  the  presence  of  glide 
and  defines  it  as  occurring  in  a  random  direction,  changing  in  a  random  manner  with  an 
initial  horizontal  velocity  of  2  to  4  feet  per  second  (fps)  (Natick  1996:4).  Under  Natick 
advisory  correspondence,  the  modeling  of  the  glide  behavior  is  made  to  mimic  (as  best  as 
possible)  the  true  behavior  of  glide  under  similar  conditions.  When  the  paratrooper 
object  has  reached  steady  state  descent,  the  onset  of  glide  is  modeled  fi-om  a  uniform 
draw  between  0  and  360  degrees  for  direction  and  between  0  to  4  fps  for  velocity. 
Subsequent  changes  in  the  glide  follows  the  initial  direction  with  little  deviation,  while 
glide  continues  to  vary  between  0  to  4  fps  (Watkins  1996). 

With  the  inclusion  of  random  winds,  the  Purvis  method  of  trajectory  propagation 
reveals  limitations  which  does  not  vitiate  the  model  by  any  means.  In  the  presence  of 
winds,  the  paratrooper  object  takes  on  the  direction  and  velocity  of  the  winds 
immediately.  Additionally,  in  the  presence  of  winds,  the  model  becomes  highly  sensitive 
to  the  propagation  time  step.  A  trade  off  exists  \vith  time  step  and  simulation  run  length, 
where  the  smaller  the  time  step  the  longer  a  simulation  run  takes  until  completion.  The 
section  on  verification  and  validation  discusses  further  the  effects  of  different  time  steps. 
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3.3  Verification  and  Validation 

Both  verification  and  validation  are  continual  processes  along  every  stage  of 
simulation  model  development.  Verification  involves  the  determination  that  a  simulation 
model  is  performing  as  the  developer  intends  while  validation  concerns  whether  or  not 
the  conceptual  model  on  which  the  simulation  model  is  based  accurately  represents  the 
system  being  studied  (Law  et  al.  1991:299). 

3.3.1  Verification 

A  critical  milestone  in  verification  is  the  correct  translation  of  the  Purvis 
FORTRAN  Model  into  MODSIM  III.  Due  to  the  vast  amounts  of  information  provided 
by  both  models,  a  graphical  form  of  verification  is  employed  using  MATLAB  graphing 
capabilities.  Table  4  displays  the  largest  absolute  error  in  the  state  variables  of  the 
paratrooper  object  when  compared  to  the  state  variables  of  the  FORTRAN  model.  Based 
on  these  small  errors,  we  conclude  that  the  Purvis  FORTRAN  Model  has  been  correctly 
translated  into  MODSIM  III,  since  the  errors  are  small  enough  to  be  attributable  to 
computational  differences  in  the  hardware  or  software  environment. 
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Table  4 

Error  Between  FORTRAN  and  MODSIM  III  Models 


Paratrooper  State 
Variables 

Largest  Absolute  Error 

Altitude  (ft) 

0.1827 

Down  Range  (ft) 

0.2840 

Off  Range  (ft) 

0.0000 

Velocity  (fps) 

0.3990 

Airspeed  (fps) 

0.3990 

Mach  Number 

0.0003 

Dynamic  Pressure  (psf) 

0.2540 

Axial  Acceleration  (gees) 

0.0294 

Trajectory  Angle  (deg) 

0.0143 

Pitch  Angle  (deg) 

0.4538 

Alpha  (deg 

0.4513 

Drag  Area  (ft^) 

0.0000 

In  integrating  the  paratrooper  objects  with  the  Petry  C-17  aircraft  and  vortex 
objects,  the  interaction  of  the  objects  is  another  area  where  verification  conies  into  play. 
Of  particular  interest  are  the  greenLight  and  pollVortices  Methods',  i.e.,  the  interactions  of 
the  aircraft/paratrooper  and  paratrooper/wake  vortex,  respectively.  The  individual 
paratrooper  objects  keep  track  of  all  relational  information;  e.g.,  data  on  which  aircraft  a 
particular  paratrooper  object  jumped  from,  the  time  of  exit,  the  (x,  y,  z)-coordinate  of  the 
exit  point,  and  airspeed  at  time  of  exit  are  all  maintained  within  the  paratrooper  object. 
Thus,  verification  is  greatly  facilitated.  The  paratrooper  objects  also  keep  track  of  all 
relational  information  regarding  the  wake  vortices  from  upstream  aircraft.  Again,  each 
paratrooper  object  tracks  all  known  wake  vortex  core  locations  and  computes  relational 
distances.  Verification,  although  cumbersome  when  several  paratrooper  objects  and 
wake  vortex  objects  are  in  the  same  airspace  at  the  same  time,  is  facilitated  with  the 
capability  of  the  paratrooper  objects  to  display  its  search  of,  and  distance  from,  the  wake 
vortices. 
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During  large  preliminary  test  runs  with  multiple  aircraft  and  several  hundred 
paratrooper  objects,  the  need  to  either  reduce  run  time  or  acquire  a  faster  computer 
becomes  apparent.  The  addition  of  more  aircraft  and  more  jumpers  further  reduces  the 
speed  of  the  simulation  model.  The  obvious  choice  of  focus  is  the  trajectory  propagation 
time  step,  dt,  of  the  paratrooper  objects.  As  in  the  Purvis  Model,  the  paratrooper  object’s 
original  dt  is  0.0005  of  a  second;  every  paratrooper  object  calculates  its  new  trajectory 
position  every  0.0005  seconds!  Sensitivity  analysis  done  on  different  dts  show  that 
certain  system  state  variables  of  the  paratrooper  object  are  highly  sensitive  to  the  dts 
under  different  wind  conditions.  Two  configurations  are  considered— no  winds,  and 
head/cross  winds.  In  both  cases,  increasing  the  dt  can  make  the  trajectory  propagation 
calculations  unstable,  usually  starting  with  the  rotational  angles  and  spreading  to  all  other 
state  variables. 

Larger  Propagation  Time  Step/No  Winds.  With  no  winds  present,  three 
propagation  time  steps  are  used  in  finding  the  sensitivity  of  model  stability  to  changes  in 
the  time  step  using  the  same  parachute-payload  system  configuration  as  the  original 
Purvis  Model.  The  three  dts  are  0.0005,  0.001,  and  0.05.  As  the  dt  increases,  run  length 
decreases,  although  larger  absolute  errors  become  apparent  as  the  propagation  scheme 
becomes  unstable  (particularly  in  the  rotational  state  variables).  Table  5  shows  the 
maximum  absolute  error  experienced  at  each  of  the  dts.  The  large  errors  experienced  in 
the  descent  velocity  and  airspeed  are  not  due  to  larger  values  but  rather  can  be  attributed 
to  an  earlier  inflation  of  the  canopy,  thus  causing  the  parachute-payload  system  to 
decelerate  sooner  at  the  larger  dts.  The  descent  velocity  and  the  airspeed  are  the  same 
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once  steady  state  is  reached  in  all  three  propagation  time  steps.  The  instability  in  the 
rotational  state  variables  is  readily  seen  graphically  as  dt  is  decreased  (Appendix  E).  The 
positional  state  variables,  however,  still  have  relatively  small  absolute  errors  associated 
with  decreased  dts,  with  the  final  (x,  y,  z)-coordinate  fairly  close  to  the  original  location 
using  the  original  dt  of  0.0005.  This  is  an  important  observation,  as  shorter  run  length  is 
desired. 

The  same  experimentation  is  done  on  the  parachute-payload  system  configured  as 
a  paratrooper  object.  Again,  similar  results  are  experienced  using  the  same  three  dts.  A 
favorable  run  time  length  is  gained  using  the  larger  propagation  time  step  of  0.05 
seconds,  and  hence  becomes  a  candidate  for  use  in  the  paratrooper  objects  with  the 
knowledge  of  the  instability  experienced  by  the  rotational  state  variables  and  an  earlier 
canopy  inflation.  The  selection  of  this  dt  is  conditional  on  system  performance  under  the 
case  where  head/cross  winds  are  present. 

Larger  Propagation  Time  StepAVinds.  Further  experiments  under  different  dts 
are  made  under  the  conditions  with  winds  present,  again  using  the  same  parachute- 
payload  system  as  in  the  original  Purvis  Model.  A  constant  crosswind  velocity  of  5  fps  is 
used  to  determine  the  effects  of  the  increased  dts  on  the  system  state  variables.  Similar 
results  of  the  case  with  no  winds  are  observed;  however,  now  with  error  in  the  y-axis, 
increased  instability  is  experienced  in  the  rotational  state  variables.  Table  5  shows  the 
maximum  absolute  error  experienced  under  wind  conditions  with  the  increased  dts. 
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Table  5 

Max  Error  From  Changes  In  Propagation  Time  Step 


System 

State 

Variables 

No  Winds 

Winds 

0.0005 

dt 

0.001 

0.05 

0.0005 

dt 

0.001 

0.05 

Altitude  (ft) 

0.182 

10.166 

IffpM 

0.149 

11.203 

Down  Range  (ft) 

WBi 

0.806 

23.251 

0.811 

23.250 

Off  Range  (ft) 

mmt 

0.000 

0.000 

0.143 

0.131 

1.000 

Velocity  (fps) 

0.399 

0.332 

25.558 

0.389 

0.366 

25.542 

Airspeed  (fps) 

0.399 

0.332 

25.558 

0.389 

0.328 

25.542 

Mach  Number 

0.000 

0.000 

0.023 

0.000 

0.000 

0.023 

Dynamic  Pressure  (psf) 

0.254 

0.153 

22.731 

0.258 

0.148 

22.817 

Axial  Acceleration  (gees) 

0.029 

0.170 

38.560 

0.029 

0.171 

38.578 

Trajectory  Angle  (deg) 

0.014 

0.094 

1.310 

0.019 

0.095 

1.304 

Pitch  Angle  (deg) 

0.454 

19.995 

176.378 

0.523 

19.696 

152.811 

Alpha  (deg 

0.451 

19.981 

202.024 

0.477 

19.655 

202.737 

Drag  Area  (ft^) 

0.000 

0.000 

4.998 

0.003 

0.003 

4.998 

When  applied  to  the  paratrooper  objects  under  wind  conditions,  the  entire 
propagation  scheme  becomes  uncontrollably  unstable  early  in  the  trajectory  life  of  the 
paratrooper  object  prior  to  the  inflation  of  the  canopy.  After  a  mechanistic  breakdown  of 
the  trajectory  calculations,  it  turns  out  that  this  instability  has  always  been  present,  but  is 
muted  due  to  the  null  value  of  wind  effects  and  is  therefore  zeroed  out.  With  the  presence 
of  winds,  this  instability  manifests  earlier  in  the  trajectory  as  the  dt  increases. 

In  the  pursuit  to  achieve  a  more  efficient  run  length,  it  becomes  necessary  to  use 
two  separate  propagation  time  steps  at  different  stages  of  the  paratrooper  object’s 
trajectory.  Through  trial  and  error,  the  large  time  step  of  0.001  seconds  prior  to  full 
canopy  inflation  and  the  smaller  time  step  of  0.01  seconds  after  full  canopy  inflation 
gives  a  reasonable  solution  to  a  shorter  run  length  while  maintaining  stability. 
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3.3.2  Validation 


Under  the  guise  of  object-oriented  programming,  the  paratrooper  object  is  the 
Purvis  Model.  Henceforth,  validation  centers  on  the  configuration  from  the  parachute- 
payload  system  used  in  the  Purvis  Model  into  a  combat-gear  outfitted  paratrooper  using  a 
T-IOC.  Our  validation  first  looks  into  the  distribution  of  a  representative  sample  of 
paratrooper  weights.  Secondly,  a  graphical  comparison  between  actual  paratrooper 
trajectories  and  trajectories  generated  by  the  integrated  airdrop  model  is  accomplished. 
Finally,  Petry  (1997)  compares  actual  testing  results  from  Edwards  AFB  and  Fort  Bragg 
to  results  generated  by  the  integrated  airdrop  model. 

Paratrooper  Weights.  The  jumper  manifest  from  C-141  Jumper-Head  to  D-Bag 
Clearance  testing  is  used  as  a  representative  sample  in  finding  a  distribution  for 
paratrooper  weights.  The  paratroopers  in  this  testing  are  outfitted  with  combat  gear  and 
use  the  T-IOC  parachute  and  reserve.  Using  the  software  package  BestFit,  any  one  of  the 
distributions  in  Table  6  makes  an  appropriate  fit  using  the  Kolmogorov-Smimov  (K-S) 
Test  Statistic  with  a  level  of  significance  of  a  =  0.05,  and  n  =  82. 


Table  6 

Fitted  Distributions  and  Associated  K-S  Test  Statistics 


Distribution 

Rank 

K-S  Test  Statistic 

Logistic(247,  13.34) 

1 

0.04008 

2 

0.06968 

Beta(5.35,  4.87)  X  151  +  167 

3 

0.07252 

x\2A6) 

4 

0.07974 

The  normal  distribution  may  be  the  easiest  and  most  readily  recognizable 
distribution  from  among  the  top  choices  with  good  fits.  By  examining  both  the 
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probability-probability  (P-P)  and  the  quantile-quantile  (Q-Q)  plots,  the  use  of  the  normal 
distribution  is  a  good  choice  as  the  representative  distribution  for  paratrooper  weights. 


Figure  7.  Normal  P-P  Plot  of  Paratrooper  Weights 


Figure  8.  Normal  Q-Q  Plot  of  Paratrooper  Weights 
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Cinetheodolite  Data  Comparison.  A  method  of  recording  actual  paratrooper 
trajectories  is  the  cinetheodolite  (Cine-T)  method,  which  uses  six  synchronized  cameras 
recording  the  same  paratroopers  from  different  angles.  This  itself  is  a  crude  recording 
method  subject  to  human  error  during  the  post-processing  data  collection,  which  is 
essentially  watching  the  recorded  video  and  extrapolating  the  paratrooper  positions  from 
frame-by-frame  playback  of  the  recorded  trajectory  (Dassow  1997).  According  to 
Dassow  the  positions  are  extrapolated  using  the  known  angles  of  the  individual  cameras 
and  triangulating  using  the  information  from  the  other  cameras.  A  further  limitation  of 
the  Cine-T  is  that  trajectories  are  not  recorded  to  impact  on  the  DZ.  Using  the  Cine-T 
data  from  actual  jumps,  only  a  cursory  visual  comparison  can  be  done  between  actual 
trajectories  recorded  from  paratrooper  jumps  using  the  C-17  as  the  jump  platform,  and 
trajectories  generated  by  the  simulation  model.  This  comparison  is  essentially  the 
implementation  of  a  Turing  Test.  Although  the  exact  conditions  surrounding  the  Cine-T 
jumps  in  Figure  9  are  not  known  (i.e.,  wind  conditions  and  airspeed  of  the  aircraft)  and 
hence  carmot  be  duplicated,  valuable  information  is  still  gained  from  this  comparison. 
Judging  from  the  trajectories  of  the  Cine-T  trajectories,  a  reasonable  assumption  is  made 
that  both  head  wind  and  cross  wind  are  present,  and  subsequently  are  incorporated  in  the 
generation  of  the  trajectories  seen  in  Figure  10. 
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Furthermore,  given  the  current  assumptions  of  the  simulation,  further  refinement  in 
trajectory  generation  suggests  using  Turing  test  for  further  validation. 

Encounter  Rate  Comparisons.  Petry  (1997)  compares  actual  vortex  encounter 
rate  test  results  from  Edwards  AFB  (EAFB)  testing  with  results  generated  by  the 
simulation  model. 


Table  7 

Results  Comparisons  of  Simulation  and  Actual  EAFB  Test 


Simulation 

Flight 

Mean 

0.13500 

0.16250 

Variance 

0.03171 

0.03470 

Observations 

50.00000 

20.00000 

Pooled  Variance 

0.03255 

Hypothesized  Mean  Difference 

0.00000 

df 

68.00000 

t  Stat 

-0.57611 

P(T<=t)  one-tail 

0.28322 

t  Critical  one-tail 

1.66757 

P(T<=t)  two-tail 

0.56644 

t  Critical  two-tail 

1.99547 

3.4  Limitations 

With  the  integration  of  the  Petry  C-17  aircraft  and  vortex  objects,  the  limitations 
of  the  paratrooper  objects  define  the  limitations  of  the  simulation.  One  of  the  major 
assumptions  of  all  parachute/payload  system  models  encountered  during  the  research 
stages  of  this  thesis  effort  is  that  the  payload  is  non-intelligent  or  “unconscious.”  This  is 
also  the  case  with  the  paratrooper  objects.  In  reality,  paratroopers  train  to  control  the 
direction  of  flight  under  a  T-IOC  canopy.  This  is  an  important  limitation  in  that  with 
non-intelligent  paratrooper  objects,  all  trajectory  and  scatter  information  only  provides  at 
best  an  upper  bound,  or  worst  case  performance  measure.  Furthermore,  starting  from 
exit,  other  limitations  are  as  follows.  First,  the  right  and  left  paratrooper  objects  exit  the 
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C-17  aircraft  object  at  discrete,  time  intervals  every  0.5  seconds.  Second,  paratrooper 
objects  immediately  assume  the  velocity  conditions  of  the  head/cross  winds.  This 
relationship  tends  to  overestimate  actual  wind  effects  on  the  paratrooper  objects,  and  thus, 
the  paratrooper  objects’  impact  points  are  again  overly  conservative.  Third,  the 
paratrooper  trajectory  propagation  scheme  in  not  affected  by  the  presence  of  wake 
vortices.  The  paratrooper  objects  do  not  interact  with  the  wake  vortices,  thus  causing  the 
impact  points  to  be  determined  without  taking  into  account  the  effect  of  an  encounter. 
Fourth,  the  modeling  of  the  paratrooper  objects  assumes  that  the  Earth  is  flat;  thus,  all 
paratrooper  objects  land  at  the  same  ground  altitude. 
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4.  ANALYSIS 


In  the  development  of  the  C-17  Paratrooper/Wake  Vortex  Simulation  Model,  the 
primary  measure  of  effectiveness  (MOE)  is  the  encounter  rate  between  paratrooper  and 
vortex  objects  (Petry  1997)  which  the  model  provides.  The  model,  however,  also 
provides  other  information  that  can  be  subjected  to  post-processing  procedures.  The 
model  provides  DZ  dispersion  information  of  the  paratroopers— information  that  is  useful 
in  the  Landing  Plan  of  the  backward  planning  sequence  used  by  airborne  commanders. 
First  and  foremost,  the  Landing  Plan  must  support  the  Ground  Tactical  Plan  and  be 
supported  by  the  Air  Movement  Plan.  Elements  of  the  Landing  Plan  are  found  in  both  the 
Ground  Tactical  Plan  and  the  Air  Movement  Plan.  Embedded  in  the  Landing  Plan,  an 
Assembly  Plan  is  generated.  Data  provided  by  the  simulation  model  on  DZ  dispersion  is 
useful  for  the  Assembly  Plan.  With  guidance  from  Klimack  (1997)  and  the  82"*^  ABN 
DIV  ASOP,  a  simplified,  though  not  trivialized,  airborne  drop  is  constructed  using 
Nijmegen  DZ  at  Fort  Bragg.  Two  MOEs  are  used:  (1)  mean  paratrooper  distance  from 
assembly  area  (AA);  and,  (2)  distribution  of  the  dispersion  across  the  width  as  well  as 
down  the  length  of  the  DZ. 

4.1  The  Scenario 

Nijmegen  DZ  is  one  of  the  smaller  DZs  used  at  Fort  Bragg.  It  measures  4950  feet 
long  and  3000  feet  wide.  Nijmegen  DZ  has  a  personnel  point  of  impact  (PI)  at  1050  feet 
into  the  DZ,  meaning  green  light  is  not  lit  until  there  is  high  confidence  that  paratroopers 
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will  land  at  least  1050  feet  into  the  DZ.  When  a  C-17  flies  at  135  knots,  Nijmegen  is 
referred  to  as  a  17  second  DZ.  When  paratroopers  exit  the  aircraft  at  0.5-second 
intervals,  this  allows  for  a  stick  size  of  1 8  jumpers  with  one  stick  exiting  each  door  per 
pass. 

Three  airborne  rifle  companies  are  used  with  three  separate  AAs  for  each 
company.  An  AA  for  the  airborne  unit  is  the  place  where,  upon  landing,  the  unit  is 
tactically  organized  and  ready  to  fight  (82"‘'  ABN  DIV  ASOP  1985:3-35  -  3-45).  The  AA 
should  be  as  close  as  possible  to  where  the  paratroopers  will  land  in  the  DZ,  with  the 
ASOP  suggesting  areas  that  lie  along  the  flanks  (instead  of  the  ends)  of  the  DZ.  The  three 
AAs  chosen  for  this  scenario  all  lie  on  the  very  edge  of  the  DZ.  The  first  AA  is  to  port  of 
the  flight  path  at  the  PI.  The  second  AA  is  located  2000  feet  down  range  from  the  PI  and 
to  starboard  of  the  flight  path.  The  last  AA  is  located  4000  feet  down  range  from  the  PI, 
and  is  also  to  starboard  of  the  flight  path.  Figure  1 1  depicts  the  scenario  used. 

With  three  companies  assembling  at  three  different  areas,  the  concept  of  cross 
loading  is  essential  to  the  successful  execution  of  the  Assembly  Plan.  When  the  C-17  is 
loaded,  the  different  units  are  partitioned  and  loaded  among  the  aircraft  such  that  the 
paratroopers  land  according  to  their  designated  AAs.  Cross-loading  both  enhances  unit 
survivability,  and  maintains  the  tactical  integrity  of  the  operation.  Cross  loading  can  also 
expedite  assembly  if  done  correctly.  It  is  in  this  Assembly  Plan  where  the  mean  distances 
from  the  AA  and  the  dispersal  distributions  on  the  DZ  are  used  as  MOEs. 
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4.2  The  Design 

Similar  to  the  Petry  experimental  design  done  with  encounter  rates,  two  features 


of  aircraft  formation  geometry  are  used  in  finding  the  MOEs:  trail  distance  and  lateral 


distance.  In  other  words,  how  do  these  two  features  affect  the  MOEs? 


Table  8 

Design  Point  Description 


Design  Point 

Seed  Used 

Trail  Distance  (ft) 

Lateral  Separation  (ft) 

Real 

Coded 

Real 

Coded 

1 

1 

- 

- 

2 

2 

BBbhiW 

- 

500 

+ 

3 

3 

32,000 

+ 

0 

- 

4 

4 

32,000 

+ 

500 

+ 

5 

5 

23,500 

0 

250 

0 

6 

6 

15,000 

- 

250 

0 

7 

7 

23,500 

0 

0 

- 

8 

8 

32,000 

+ 

250 

0 

9 

9 

23,500 

0 

500 

+ 
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In  approaching  this  experimental  design,  multiple  replications  of  the  model  are  performed 
at  different  set  values  for  the  trail  and  lateral  distances.  Table  8  defines  the  design  points 
used  while  Figure  12  shows  how  both  the  trail  and  lateral  distances  are  varied.  A  simple 
2^  factorial  design  is  used  with  center  point  runs.  However,  additional  points  of  interest 
are  also  considered.  The  additional  face  points  are  included  in  the  design  to  gain  better 
insight  while  staying  within  our  region  of  operation.  The  idea  of  adding  axial  points  with 
ad>  1.0  extends  the  region  in  areas  which  would  not  provide  any  added  benefit  {ad  is  the 
distance,  in  coded  terms,  from  the  center  point  and  is  different  from  the  a  in  statistical 
significance  levels).  An  axial  point  for  lateral  spacing,  for  example,  merely  places  the 
trail  aircraft  on  the  opposite  side  of  the  lead  aircraft,  an  area  that  is  already  accounted  for 


Figure  12.  Variations  On  Formation  Geometry 
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by  center-point  observations.  In  short,  the  design  is  a  central  composite  structure  with  ad 
=  1.0,  or  a  face-centered  design  (Neter  et  al.  1996:1282-1286). 
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Figure  13.  Face-Centered  Design 

In  the  design,  the  distance  of  each  individual  paratrooper  impact  point  to  its 
respective  AA  is  not  looked  at  as  an  individual  observation.  The  mean  distance  from  the 
AAs  of  all  the  paratroopers  in  the  entire  two-ship  formation  defines  the  MOE  for  a  single 
run.  It  is  similar  to  the  principle  applied  when  looking  at  the  distribution  on  the  DZ.  An 
individual  impact  point  carmot  give  meaningful  information  on  the  distribution.  It  is  only 
when  the  dispersion  is  looked  at  in  aggregate  that  a  meaningful  distribution  can  be  found. 
(This  is  similar  to  the  method  of  batch  means  used  for  non-terminating  simulations, 
where  the  batch  size  is  the  number  of  paratroopers  used  in  each  of  the  replications.)  Each 
run  has  an  associated  mean  distance  per  company  T,  and  sample  variance  5/  where 
Yim  =  distance  from  designated  AA  for  company  i 
c  =  number  of  runs  per  design  point 
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i 

J 

k 

n 

and 


=  designator  for  company  number  and  respective  AA 
=  run  number 
=  the  design  point 
=  number  of  paratroopers 


and 


n  —  1  m=l 

Furthermore,  each  design  point  has  an  associated  overall  mean  distance  7*,  such  that 

n=-tYy 

C  jai 

with  variance 

s- 

C  —  1  y  =  I 

—2 

and  mean  variance  Ski  where 


which  itself  has  variance 


1 


S-2  -  .  ^ 

C  -  1  y=| 

The  dispersion  of  the  paratroopers  across  the  width  and  dovra  the  length  of  the  DZ 
when  fitted  to  a  distribution  as  a  paired  set  defines  the  second  MOE  of  interest  to  the 
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airborne  commanders  (Klimack  1997).  Although  the  82"'*  ABN  DIV  ASOP  does  not 
specifically  define  a  favored  distribution  for  these  dispersions,  it  is  beneficial  to  see 
whether  the  dispersion  on  the  DZ  supports  the  Ground  Tactical  Plan.  It  is  desirable  that 
the  paratroopers  disperse  around  the  flight  path  of  the  aircraft  with  some  form  of 
distribution.  Intuitively,  the  dispersion  across  the  width  of  the  DZ  may  depend  greatly  on 
the  lateral  spacing  between  the  aircraft.  These  paired  MOEs  are  extracted  during  post¬ 
processing  of  the  design  results. 

The  initial  number  of  replications  used  at  each  design  point  is  twelve.  In 
considering  the  number  of  runs  needed,  the  mean  distances  from  each  of  the  three  AAs 
are  used  with  a  very  conservative  tolerance  of  5  feet.  If  the  need  to  create  more 
replications  becomes  evident,  performing  further  runs  is  a  simple  task  to  do.  (Post¬ 
processing  reveals  that  the  12  replications  are  more  than  enough  to  be  within  the  5  feet 
tolerance  at  an  alpha  level  of  a=0. 05.) 

4.3  The  Results 

4.3.1  Mean  Distance  From  Assembly  Areas 

Table  9  summarizes  the  results  of  the  12  runs  at  each  design  points.  It  contains 

only  the  Yki  and  sli  of  each  design  point  for  each  of  the  three  companies.  Table  24  in 

Appendix  F  shows  the  complete  results  with  standard  deviations  for  each  of  the  Yki  and 
—2 

Ski ,  including  fairly  small  variations  in  all  the  variables.  However,  looking  at  the  results 
within  each  design  point,  the  variation  in  mean  distance  is  consistent  with  our  intuition— 
that  lateral  spacing  has  a  large  influence  while  trail  distance  has  a  small  effect  on  mean 
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Table  9 

Summary  Of  Results  For 
Mean  Distance*  To  Assembly  Areas 


Design  Poinf 

1 

Company 

2 

3 

1 

436.228 

362.8717 

387.8438 

74.0577 

34.5662 

57.795 

2 

504.9642 

287.2827 

319.209 

(-,+) 

98.2338 

83.8582 

96.8437 

3 

428.0379 

361.4047 

392.7657 

(+,-) 

72.9951 

30.6491 

61.6572 

4 

497.1976 

285.635 

327.4557 

(+,+) 

96.2339 

84.306 

101.5088 

5 

468.5109 

319.6805 

351.9928 

(0,0) 

77.5729 

51.1272 

72.7923 

6 

468.2581 

323.4957 

352.2123 

(-,0) 

79.2695 

52.4096 

73.2678 

7 

433.0555 

359.7963 

391.7037 

(0,-) 

74.6482 

31.3088 

58.2386 

8 

460.9896 

322.2542 

360.2295 

(+,0) 

77.7871 

50.4006 

76.3632 

9 

499.979 

286.9205 

322.376 

(0,+) 

97.8993 

85.2177 

100.2291 

'  Denotes  coded  factor  levels  at  each  design  point. 

*  Mean  distance  is  above  its  standard  deviation.  Standard  deviation  is  italicized. 


distance.  This  observation  is  supported  by  looking  at  the  correlation  matrix,  Table  10,  of 


the  trail  {Xj)  and  lateral  {Xi)  distances  and  all  the  7,  ’s. 


Table  10 

Correlation  Matrix 


Variables 

-^2 

Yx 

Yi 

F3 

^1 

1.0000 

^2 

0.0000 

1.0000 

Ki 

-0.1112 

0.9811 

1.0000 

Yi 

-0.0192 

-0.9887 

-0.9823 

1.0000 

Yi 

0.1023 

-0.9811 

-0.9895 

0.9839 

1.0000 

A  more  sophisticated  method  of  looking  at  the  effects  of  changes  in  trail  and 
lateral  distance  can  be  accomplished  through  a  RSM  approach.  With  the  central 
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Table  11 

Response  Surface  Parameter  Estimates 


Yx 

Y2 

Yi 

Intercept 

466.7435 

320.6825 

353.9703 

X, 

-3.8709 

-0.7260 

3.5310 

X2 

34.1366 

-37.3725 

-33.8788 

x,^ 

-1.2359 

1.6915 

1.2618 

X,X2 

0.1059 

-0.0452 

0.8312 

0.6575 

2.1750 

2.0808 

composite  design  used,  the  comer  points  provide  estimations  for  the  linear  main  effects 


and  interaction  effects,  the  axial  points  provide  an  estimation  for  the  quadratic  main 


effects,  and  the  center  points  provides  both  a  pure  error  estimate  for  and  allows  for  a 


lack  of  fit  test  (Neter  et  al.  1996:1282).  The  full  functional  form  of  the  response  surface 


that  is  generated  with  coded  variables  is  provided  in  Table  11  with  the  analysis  of 


variance  (ANOVA)  provided  in  Table  12. 


Statistical  significance  in  the  parameter  estimates  can  be  used  in  finding  a 


parsimonious  model;  however,  with  an  =  1.0  tests  of  significance  for  the  quadratic 


terms  cannot  be  performed  (Neter  et  al.  1996:  1286).  Therefore,  all  the  parameters  are 


kept.  For  Y\ ,  the  intercept,  X,,  and  X2  are  statistically  significant.  For  Yi ,  only  the 


intercept  and  are  significant.  Finally,  T3  has  the  same  significant  effects  as  Y\ .  In  all 


Table  12 

Analysis  of  Variance  for  Significance  of  Regression  In  Multiple  Regression 


Response 

Source 

DF 

Sum  of  Squares 

Mean  Square 

F-Ratio 

Prob  >  F 

Model 

5 

85028.350 

17005.700 

808.795 

<0.0001 

h 

Error 

2144.644 

21.000 

Total 

Bl 

87172.995 

n 

100782.640 

20156.500 

983.682 

<0.0001 

h 

mm 

2090.07 

20.500 

Total 

WSM 

102872.710 

Model 

5 

83713.418 

16742.500 

800.166 

<0.0001 

Yi 

Error 

102 

2134.225 

20.900 

Total 

107 

85846.642 

three  cases,  values  >  0.97  are  experienced.  The  graphical  representations  of  these 
functions  follow. 


LATERAL  DISTANCE 


TRAIL  DISTANCE 


Figure  14.  Response  Surface  Plot  for  Mean  Distance 
From  AA  1 


TRAIL  DISTANCE 


LATERAL  DISTANCE 


Figure  15.  Response  Surface  for  Mean  Distance 
From  AA  2 


TRAIL  DISTANCE 


LATERAL  DISTANCE 


Figure  16.  Response  Surface  for  Mean  Distance 
From  A  A  3 


The  surface  plots  of  the  mean  distances  for  each  of  the  AAs  clearly  show  the  small  effect 
which  changes  in  trail  distance  have  compared  to  the  effects  of  changes  in  lateral 
distance. 


4.3.2  Distribution  On  The  Drop  Zone 

In  finding  a  representative  distribution  of  paratrooper  dispersal  across  the  width 
and  down  the  length  of  a  DZ,  the  primary  focus  of  data  description  can  easily  be  lost 
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within  the  post-processing  investigative  analysis.  This  occurs  during  the  repetitive  nature 
of  post-processing,  as  there  are  nine  design  points  to  be  examined,  with  each  design  point 
having  12  replications  and  each  replication  having  two  variables  to  work  with. 
Representing  the  MOEs  with  2 1 6  separate  distributions  can  distort  the  underlying  reason 
for  the  investigation  in  the  first  place — ^to  find  a  descriptive  distribution  for  dispersion. 

Much  information  can  be  gathered  by  taking  all  the  impact  locations  of  the 
paratroopers  used  in  each  design  point  and  analyzing  them  as  a  whole.  Even  a  cursory 
look  at  the  basic  statistics  of  the  impact  locations  in  Table  13  yields  valuable  insight  into 
paratrooper  performance.  In  Table  13,  location  down  the  length  of  the  DZ  is  represented 
by  the  x-component  while  location  across  the  width  is  represented  by  the  y-component. 

Several  observations  stand  out.  First,  all  the  summary  statistics  for  the  x- 
components  at  each  design  point  are  similar.  Secondly,  the  means  of  the  y-components 
are  grouped  around  the  midpoint  between  the  lead  and  trail  aircraft  for  all  three  lateral 
separations  used.  Finally,  as  lateral  separation  increase,  the  variance  of  the  y-  component 
also  increases. 
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Table  13 

Summary  Statistics  Of  Dispersion  Data 


Design  Point 

Component 

Mean 

Standard  Deviation 

Min 

Max 

1 

X 

3270.4 

1191.1 

1170.9 

5371.3 

(-,-) 

y 

0.9 

59.0 

-116.9 

127.1 

2 

X 

3270.3 

1191.3 

1174.9 

5377.7 

(-,+) 

y 

251.0 

254.3 

-123.8 

619.3 

3 

X 

3272.7 

1193.6 

1090.4 

5327.3 

(+,-) 

y 

-0.2 

58.2 

-ni.i 

130.1 

4 

X 

3223.5 

1193.4 

1127.7 

5340.3 

(+,+) 

y 

250.4 

255.4 

-114.6 

612.9 

5 

X 

3247.8 

1191.2 

1128.9 

5352.8 

(0,0) 

y 

125.1 

135.7 

-113.3 

367.1 

6 

X 

3272.0 

1193.0 

1161.7 

5369.6 

(-,0) 

y 

123.4 

137.0 

-128.8 

364.2 

7 

X 

3245.8 

1190.3 

1116.4 

5349.7 

(0,-) 

y 

1.3 

57.0 

-125.3 

118.3 

8 

X 

3221.4 

1190.8 

1125.1 

5333.0 

(+,0) 

y 

123.1 

138.5 

-119.2 

365.7 

9 

X 

3247.1 

1194.2 

1129.6 

5348.8 

(0,+) 

_ Z _ 

249.2 

258.8 

-110.4 

614.4 

Interpreting  these  statistics  in  relation  to  the  first  MOE  shows  that  the  two  results 
support  one  another.  For  the  first  and  second  MOE,  trail  distance  has  little  effect  whereas 
lateral  distance  influences  the  MOEs  considerably.  Another  observation  is  that  variation 
in  the  y'-component  increases  as  lateral  distance  increases,  which  occurs  because  the 
lateral  separation  scatters  the  paratroopers  in  a  wider  path  down  the  length  of  the  DZ.  In 
all  cases,  each  design  point  generates  dispersion  down  the  length  of  the  DZ  that  appears 
to  be  uniformly  distributed  regardless  of  trail  or  lateral  separation  (Appendix  G). 
However,  dispersion  across  the  DZ  is  highly  dependent  on  the  lateral  separation.  With  no 
lateral  separation,  the  dispersal  distribution  is  heavily  skewed  towards  either  side  of  the 
flight  path.  As  lateral  separation  increases  a  bi-modal  effect  occurs,  creating  two 
distinguishable  distributions  each  exhibiting  characteristics  similar  to  the  distribution  of 
the  no  lateral  separation  scenario.  Figures  17,  18,  and  19  demonstrate  this  bi-modal 
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ProporUon 


distribution  and  the  separation  effect  experienced  with  lateral  separation  at  all  three  trail 
distances. 


Figure  17.  Dispersion  Distribution  Across  DZ  With  No  Lateral  Separation 

I  I  5000  Fn  T rail  Distance 


Figure  18.  Dispersion  Distribution  Across  DZ  With  Lateral  Separation  of  250  Feet 
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Figure  19.  Dispersion  Distribution  Across  DZ  With  Lateral  Separation  of  500  Feet 


Although  the  airborne  commander  does  not  use  mean  distance  to  AAs  or  DZ 
dispersal  distribution  directly,  the  commander  wants  to  deliver  the  airborne  units  onto  the 
DZ  in  the  best  configuration  possible,  which  is,  in  part,  a  function  of  the  mean  distance 
and  the  distribution  of  the  paratrooper  dispersal  (Klimack  1997).  How  this  may  affect 
operational  planning  relates  to  the  selection  of  the  AAs  as  to  minimize  the  mean  distance 
from  the  chosen  AAs.  In  turn,  mean  distance  to  AAs  and  DZ  dispersion  are  both 
functions,  in  part,  of  aircraft  formation  geometry.  This  does  not  suggest  changes  to  the 
ASOP.  It  does,  however,  provide  a  tool  with  which  commanders  can  consider  both 
paratrooper  safety  regarding  encounter  rates  with  wake  vortices  and  DZ  assembly  in  using 
specific  operational  formations  from  which  to  jump. 
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5.  CONCLUSIONS 


5.1  Overview 

The  utility  of  object-oriented  programming  is  demonstrated  with  the  modeling  of 
the  complex  airborne  operations  environment.  Taking  a  mechanistic  view  of  the 
activities  involved  in  airborne  operations,  a  more  penetrating  look  into  a  systems  aspect  of 
relational  interactions  between  entities  involved  in  the  airborne  activities  is  possible. 
This  Thesis  effort  with  the  parallel  efforts  of  Petry  (1997)  has  done  just  that.  In 
translating  an  established  6-DOF  trajectory  model  into  an  object  and  transforming  that 
object  to  reflect  the  aerodynamic  characteristics  of  a  combat  equipped  paratrooper,  the 
mechanistic  breakdown  is  completed  at  the  level  of  resolution  desired.  In  approaching 
the  systems  view,  these  paratrooper  objects  are  integrated  with  the  Petry  C-17  aircraft  and 
vortex  objects  into  a  simulation  model  that  allows  for  a  natural  hierarchy  in  a  system 
build-up. 

5.2  Conclusions 

In  applying  the  simulation  model  to  a  simplified  airborne  operation,  three  airborne 
rifle  companies  are  used  to  investigate  the  effects  that  trail  and  lateral  distances  have  on 
(1)  mean  distance  to  AAs  and  (2)  dispersion  distribution  on  the  DZ.  In  both  cases,  trail 
distances  have  little  influence  on  the  MOEs  whereas  lateral  distances  directly  affect  them. 
Depending  on  where  the  AAs  are  chosen,  mean  distance  to  the  AA  decreases  as  lateral 
separation  increases  if  the  AA  is  chosen  on  the  same  side  of  where  separation  is  towards. 
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At  the  same  time,  it  increases  the  variance  of  the  mean  distance  because  lateral  separation 
scatters  the  paratroops  in  a  wider  path  across  the  width  of  the  DZ.  These  results  are  in 
light  of  the  one  major  assumption  that  the  paratrooper  object  is  a  non-intelligent  payload. 
For  the  airborne  commander,  this  translates  into  the  operational  parameters  of  aircraft 
formation  geometries  when  also  considering  encounter  rates  between  the  paratroopers  and 
wake  vortices  of  upstream  aircraft.  This  balancing  act  of  finding  an  operational  aircraft 
formation  which  minimizes  both  the  encounter  rate  and  the  mean  distance  to  AAs  has  a 
direct  influence  on  the  type  of  aircraft  formation  to  implement  in  addition  to  all  other 
factors  which  the  airborne  commander  considers. 

5.3  Recommendations  and  Future  Development 

The  limitations  of  this  model  can  be  used  as  the  starting  point  for  future 
development  and  improvements  on  the  objects  used.  A  distribution  can  be  found  for 
inter-exit  times  for  the  paratrooper  objects.  Along  with  this,  however,  is  the  change  in  the 
update  method  of  the  Petty  aircraft  object  for  its  positioning.  A  time  delayed  wind  effect 
can  be  incorporated  to  reduce  the  overly  conservative  influence  on  trajectory  which 
head/cross  winds  have.  Lastly,  an  object  that  defines  ground  terrain  and  elevation  can  be 
incorporated  and  made  to  interact  with  the  paratrooper  objects  to  enhance  resolution  by 
taking  into  account  natural  land  formations  such  as  DZs  with  tree  lines.  The  limitation  on 
the  dynamics  of  paratrooper/wake  vortex  interactions  in  terms  of  how  they  occur  and  the 
resulting  effects  is  best  left  for  aerospace  engineers.  However,  when  the  resources 
become  available  to  feasibly  study  the  dynamics  of  paratrooper/wake  vortex  interaction. 
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an  area  of  interest  which  needs  to  be  investigated  is  the  effects  the  interactions  may  have 
on  the  variation  in  the  paratrooper  dispersal  distributions  on  the  DZ  and  then  interpret 
these  findings  into  operational  parameters  which  can  be  easily  applied  by  the  airborne 
commander. 

One  aspect  of  object-oriented  simulation  that  can  be  applicable  to  future 
developments  is  the  idea  of  inheritance.  These  paratrooper  objects  are  modeled 
specifically  with  the  T-IOC.  An  upgrade  to  the  T-IOC  is  being  developed,  called  the 
Advanced  Tactical  Parachute  System  (ATPS),  which  may  cause  the  re-testing  of 
formation  geometries  to  investigate  whether  encounter  rates  change  with  the  change  of 
parachute  performance.  With  inheritance,  a  new  paratrooper  object  can  be  developed, 
inheriting  everything  from  the  old  paratrooper  object  but  with  the  aerodynamic 
characteristics  and  performance  of  ATPS. 

Creating  a  graphical  user  interface  (GUI)  for  the  simulation  model  can  improve 
the  model’s  utility.  A  GUI  can  also  facilitate  the  understanding  of  the  dynamics  involved 
in  the  modeling  of  the  object  interactions.  However,  if  rotational  information  is  needed, 
the  original  propagation  time  step  problem  must  be  revisited.  This  addition  to  the 
simulation  model  is  already  being  considered  with  independent  efforts  being  pursued 
outside  of  this  Thesis. 

5.4  Summary 

Without  trivializing  the  complexities  with  not  just  the  aerodynamics  involved 
with  modeling  a  parachute/payload  system  but  also  with  the  interactions  of  the 
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parachute/payload  system  with,  its  environment  as  in  airdrop  operations,  a  simplistic 
rendition  of  this  airborne  world  is  achieved  using  object-oriented  simulation.  With  the 
original  purpose  of  this  Thesis  effort  having  been  met  (to  provide  the  C-17  test  and 
evaluation  community  with  the  capability  to  assess  paratrooper  performance  during  C-17 
airdrop  formations),  it  is  applied  to  a  simplified  scenario  that  parallels  real  world  airborne 
operations  to  demonstrate  the  simulation  model’s  capabilities  and  applicability.  When 
combined  with  the  findings  of  the  Petiy  investigation  into  the  paratrooper/wake  vortex 
encounter  rates,  the  airborne  commander  can  have  better  insight  into  the  environment 
surrounding  the  operations  of  the  airborne  units. 
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APPENDIXA 

Purvis  FORTRAN  6-DOF  Code 


c 


C . SIX  DEGREE  OF  FREEDOM  FLIGHT  SIMULATION 

C .  DIRECTION  COSINE-EULER  AXES  METHOD 

C 

C  SINGLE  BODY  -  CHUTE  DECELERATION  PROGRAM 

C 

C  FOREBODY  DRAG  VS  MACH  NO.  VERSION 

C 

C . EXECUTIVE  ROUTINE 

C 

C . SOURCE  :  DR.  J.  PURVIS 

C  CCG-COURSE  F12.01 

C  MUENCHEN,  1987 

C 

C  MODIFIED  BY  DR.  K.-F.  DOHERR 

C  FOR  IBM-COMPATIBLE  PC 

C  IBM  PROFESSIONAL  FORTRAN  WAS  USED  FOR  COMPILATION 

C  8087  MATH-COPROCESSOR  IS  NEEDED 

C 

C  INPUT  FROM  DIKFILE  SIXD.DAT  EXPECTED 

C  OUTPUT  ON  DISKFILE  SIXD.GRF  AND  ON  UNIT  6  =  SCREEN 

C . 

REAL  MASS, IN, JN, MACH 
CHARACTER*!  DUM 


COMMON/INERT/MASS,XCG,XBOD, IN (3, 3) , JN (3, 3) 
COMMON/AEROS/CBAR, SAREA,MACH 

COMMON/SUBS/  XE { 3 ) , UE {3 ) , WB ( 3 ) , B ( 3, 3 ) , VWIND ( 3 ) 
COMMON/BURN/  T 

COMMON/CHUTE/  IPTS, PCDS (50) , PT (50) 

COMMON/ FBDRAG/  MPTS, PCDF(50)  ,  PM(50) 

COMMON/AEROCF/  CNA, CNA2 , CNQ, CYB, CYB2, CYR, CAO, CAA2 
COMMON/AEROCM/  CLO , CLP, CMO , CMA, CMA2 , CMQ, CNB, CNB2 , CNR 
COMMON/PROG/  DT, DTPR, TEND, HMIN 
COMMON/REF/  H, SOUND, RHO, RHOZ 

201  FORMAT(F10.3) 

202  FORMAT (Al) 

203  FORMAT (29X, 12) 

204  FORMAT(6X,F10.2,2X, F10.2) 

205  FORMAT (27X,F10. 2) 


C 

C . OPEN  FILES  FOR  INPUT  AND  OUTPUT 

C 

C . UNIT  11  =  INPUT  FILE  =  SIXD.DAT 

C . UNIT  6  =  OUTPUT  FILE  =  CON  TABULATED  OUTPUT 

C . UNIT  12  =  OUTPUT  FILE  =  SIXD.GRF  ASCII  FILE 

C 

OPEN  (UNIT  =  11,  FILE  =  'SIXD.DAT',  STATUS  =  'OLD') 
OPEN  (UNIT  =  12,  FILE  =  'SIXD.GRF') 


INPUT  SYSTEM  INERTIAL  PROPERTIES 
CALL  HEAD 

READ(11,201)  WEIGHT 
MASS=WEIGHT/32.17 
READ(11,201)  XCG 
READ(11,201)  XBOD 
DO  19  J=l,3 
DO  19  L=l,3 
JN(J,L)=0. 

IN(J,L)=0. 

19  CONTINUE 

READ(11,201)  IN(1,1) 

READ(11,201)  IN(2,2) 

READ(11,201)  IN(3,3) 

DO  2  1=1,3 
JN(I,I)=1./IN(I,I) 

2  CONTINUE 
C 


00000100 

00000200 

00000300 

00000400 

00000500 

00000600 

00000700 

00000800 

00000900 

00001000 

00001010 

00001020 

00001030 

00001040 

00001050 

00001060 

00001070 

00001080 

00001090 

00001091 

00001092 

00001093 

00001100 

00001200 

00001300 

00001400 

00001500 

00001600 

00001700 

00001800 

00001900 

00002000 

00002100 

00002200 

00002300 

00002400 

00002500 

00002600 

00002700 

00002800 

00002900 

00003000 

00003001 

00003002 

00003003 

00003010 

00003012 

00003014 

00003015 

00003020 

00003030 

00003100 

00003200 

00003300 

00003400 

00003500 

00003600 

00003700 

00003800 

00003900 

00004000 

00004100 

00004200 

00004300 

00004400 

00004500 

00004600 

00004700 
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c. 

....INPUT  INITIAL  CONDITIONS 

00004800 

c 

00004900 

DO  22  I=T,  3 

00005000 

XE(I)=0. 

00005100 

aE(I)=0. 

00005200 

WB(I)=0. 

00005300 

VWIND(I)=0. 

00005400 

22  CONTINUE 

00005500 

CALL  HEAD 

00005600 

READ(11,201)  ALT 

00005700 

READ(11,201)  HMIN 

00005800 

READ(11,201)  UE(1) 

00005900 

READ(11,201)  DE{3) 

00006000 

READ(11,201)  THETA 

00006100 

READ(11,201)  VWIND(l) 

00006200 

READ(11,201)  VWIND(2) 

00006300 

READ{11,201)  DENS 

00006400 

RHOZ=0. 002378 

00006500 

IF(DENS.EQ.O. )  GOTO  5 

00006600 

RHOZ=DENS*EXP (ALT/23111 . -.295*SIN (ALT/28860 . ) - . 213* 

00006700 

$  SIN(ALT/86580. ) ) 

00006800 

5  CONTINUE 

00006900 

c 

00007000 

c. 

.  . . . INPUT  PROGRAM  CONSTANTS 

00007100 

c 

00007200 

CALL  HEAD 

00007300 

READ(11,201)  DT 

00007400 

READ(11,201)  DTPR 

00007500 

READ(11,201)  TEND 

00007600 

c 

00007700 

c. . 

.  . . . INPUT  FOREBODY  AERODYNAMIC 

COEFFICIENTS 

00007800 

c 

00007900 

CALL  HEAD 

00008000 

READ(11,201)  CBAR 

00008100 

READ(11,201)  SAREA 

00008200 

READ(11,201)  CNA 

00008300 

READ(11,201)  CYB 

00008400 

READ(11,201)  CAA2 

00008500 

READ(11,201)  CLO 

00008600 

READ{11,201)  CLP 

00008700 

READ(11,201)  CMA 

00008800 

READ(11,201)  CMQ 

00008900 

READ(11,201)  CNB 

00009000 

READ{11,201)  CNR 

00009100 

c 

00009200 

c. . 

. . . INPUT  FOREBODY  DRAG  VS  MACH 

NO.  TABLE 

00009300 

c 

00009400 

CALL  HEAD 

00009500 

READ  (11,203)  MPTS 

00009600 

READ  (11,202)  DOM 

00009700 

DO  1  1=1, MPTS 

00009800 

READ  (11,204)  PM(I) , PCDF(I) 

00009900 

1  CONTINUE 

00010000 

c 

00010100 

c. . 

. . . INPUT  PARACHUTE  DRAG-AREA  VS  TIME  TABLE 

00010200 

c 

00010300 

CALL  HEAD 

00010400 

READ  (11,205)  DEPTIME 

00010500 

READ  (11,203)  IPTS 

00010600 

11=1 

00010700 

IF (DEPTIME. LT. 0. )  DEPTIME=0 

00010800 

IF(DEPTIME.LE.0 . )  GOTO  4 

00010900 

11=3 

00011000 

IPTS=IPTS+2 

00011100 

PT(1)=0. 

00011200 

PT(2)=DEPTIME 

00011300 

PCDS (1)=0. 

00011400 

PCDS(2)=0. 

00011500 

4  READ  (11,202)  DUM 

00011600 

DO  3  1=11, IPTS 

00011700 

READ  (11,204)  PT(I) ,PCDS(I) 

00011800 

PT(I)=PT(I)+DEPTIME 

00011900 

3  CONTINUE 

00012000 
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_  -  DO  999  I  =  1,4 

~  WRITE  (*,*)  PODS (I) 

999  CONTINUE 
C 

C . CONVERT  EULER  ANGLES  TO  DERECTION  COSINES 

C 

PSI=0. 

PHI=0 . 

RAD=l./57.295 

ST=SIN(THETA*RAD) 

CT=COS (THETA*RAD) 

SP=SIN(PSI*RAD) 

CP=COS (PSI*RAD) 

SPHI=SIN(PHI*RAD) 

CPHI=COS (PHI*RAD) 

XE(3)=-ALT 
B(l, 1)=CP*CT 
B(1,2)=SP*CT 
B (1, 3) =-ST 

B  (2, 1) =-SP*CPHI+CP*ST*SPHI 
B  (2, 2)=CP*CPHH-SP*ST*SPHI 
B (2, 3) =CT*SPHI 
B (3, 1) =SP*SPHI+CP*ST*CPHI 
B (3, 2) =-CP*SPHI+SP*ST*CPHI 
B(3,3)=CT*CPHI 
C 

C . TRAJECTORY  SIMULATION 

C 

CALL  TRAJEC 
C 

CLOSE  (11) 

CLOSE  (12) 

C 

STOP 

END 

SUBROUTINE  HEAD 
CHARACTER*!  DUM 
DO  1  1=1,5 
READ  (11,202)  DUM 
202  FORMAT (Al) 

1  CONTINUE 
RETURN 
END 

SUBROUTINE  TRAJEC 
C 

C . EQUATIONS  OF  MOTION  -  TRAJECTORY  SIMULATION 

C 

REAL  MASS, IN,  JN, MACH, MB 
INTEGER  E 

COMMON/INERT/MASS,XCG,XBOD,  IN(3,3),JN(3,3) 

COMMON/AEROS/  CBAR, SAREA, MACH 

COMMON/SUBS/  XE(3) ,UE(3) , WB (3 ) , B (3, 3) ,  VWIND ( 3) 

COMMON/BURN/  T 

COMMON /AEROF/  CN, CY, CA, CSL, CM, CSN, CDS 
COMMON/AEROV/  VPO, SALP, CALP, SBET, CBET, PB, QB,  RB 
COMMON/PROG/  DT, DTPR, TEND, HMIN 
COMMON/REF/  H, SOUND, RHO, RHOZ 
DIMENSION  BN (3, 3) , TEMP ( 3) , E (5) , DEL (3,  3) 

DIMENSION  FB(3) ,MB (3) , FE(3) , UEDOT ( 3) , WBDOT (3)  ,  BDOT(3, 3)  ,HB (3) 
DATA  E/1,2,3,  1,2/,  DEL/1  .,0.,0.,0.,1.,0.,0.,0.,1./ 

200  FORMAT (IH  , 21X, 4HDOWN, 6X, 3HOFF, 25X, 4HMACH,  3X,  7HDYNAMIC, 


00012100 

00012200 

00012300 

00012400 

00012500 

00012600 

00012700 

00012800 

00012900 

00013000 

00013100 

00013200 

00013300 

00013400 

00013500 

00013600 

00013700 

00013800 

00013900 

00014000 

00014100 

00014200 

00014300 

00014400 

00014500 

00014600 

00014610 

00014621 

00014622 

00014623 

00014700 

00014800 

00014900 

00015000 

00015100 

00015200 

00015300 

00015400 

00015500 

00015600 

00015700 

00015800 

00015900 

00016000 

00016100 

00016200 

00016300 

00016400 

00016500 

00016600 

00016700 

00016800 

00016900 

00017000 

00017100 

00017200 

00017300 

00017400 


$  3X, 5HAXIAL, 3X, lOHTRAJECTORY, 3X, 5HPITCH, 12X, 

$  4HDRAG) 

201  FORMAT (IH  , IX, 4HTIME, 3X, 8HALTITUDE, 4X, 5HRANGE, 5X,  5HRANGE, 

$  3X, 8HVELOCITY, 2X, 8HAIRSPEED, 

$  2X, 6HNUMBER, 2X, 8HPRESSORE, 2X, 

$  6HACCEL . , 4X, 5HANGLE, 6X, 5HANGLE, 3X, 5HALPHA, 4X, 

$  4HAREA) 

202  FORMATdH  ,  IX,  5H  (SEC)  ,  4X,  4H  (FT) ,  7X,  4H  (FT)  ,  6X,  4H  ( FT)  ,  5X, 

$  5H(FPS) , 5X, 5H(FPS) , 13X, 

$  5H(PSF) , 3X, 6H (GEES) , 4X,5H(DEG) , 6X, 


00017500 

00017600 

00017700 

00017800 

00017900 

00018000 

00018100 

00018200 

00018300 

00018400 
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$  5H (DEG) , 3X, 5H(DEG) , 2X, 8H (SQ. FT . ) ) 

203  FORMATdH  ,  F6 . 2 ,  IX,  5  {F8 . 1 , 2X)  ,  F6 . 2, 2X,  F8 . 1 , 2X,  F6 . 1 , 4X, 
$  F6.1,4X,F6.1,2X,F6.1,2X,F8.2) 

C 

C . FIXED  INITIAL  CONDITIONS  AND  CONSTANTS 

C 

GEES=0 . 

CDS=0.  ■ 

T=0. 

GMAX=0 . 

IEND=0 

TPR=T-DT 

TEND=TEND-0.1*DT 
ALT=-XE (3) 

H=ALT 

WRITE (12, 505)  TEND, HMIN, ALT 
WRITE (12,505)  UE(1),UE(3), VWIND ( 2 ) 

WRITE (6, 200) 

WRITE(6,201) 

WRITE (6, 202) 

C 

C . BEGIN  TRAJECTORY  LOOP 

C 

1  CONTINUE 
G=GRAV(H) 

CALL  DENS 
C 

C . COMPOTE  AERO.  VARIABLES 

C 

PB=WB(1) 

QB=WB (2) 

RB=WB(3) 

UEl=OE(l) -VWIND (1) 

OE2=UE(2)-VWIND(2) 

OE3=OE(3)-VWIND(3) 

VP=SQRT (UE1**2+UE2**2+UE3**2) 

VP0=VP 

MACH=VP/SOUND 

UB1=B(1, 1) ♦UE1+B(1,2) *UE2+B(l,3)*OE3 
UB2=B ( 2 , 1 ) *DE1+B ( 2 , 2 ) *UE2+B (2,3) *UE3 
UB3=B (3, 1) *UE1+B (3, 2) ♦UE2+B (3, 3) *UE3 
VP13=SQRT (UB1**2+UB3**2) 

C 

C . USE  SIN (ALPHA)  FOR  ALPHA  AND  SIN (BETA)  FOR  BETA 

C 

IF(VP0.LT.l.E-6)  VP0=l.E-6 

SBET=OB2/VP0 

CBET=VP13/VP0 

BETA=SBET 

IF(VP13.LT.l.E-6)  VP13=l.E-6 
SALP=UB3/VP13 
CALP=UB1/VP13 
ALPHA=SALP 
C 

C . AERODYNAMIC  AND  BODY  FORCES  AND  MOMENTS 

C 

C . ISOLATED  BODY  AERODYNAMICS 

C 

CALL  AERO 

Q=0.5*RHO*VP 

FPC=-Q*CDS 

QS= . 5*RHO*VP*VP*SAREA 
QSD=QS*CBAR 

FB ( 1 ) =-QS*CA+MASS*G*B (1, 3) +FPC*UB1 
FB (2 ) =QS*CY+MASS*G*B (2, 3) +FPC*UB2 
FB ( 3 ) =-QS*CN+MASS*G*B (3,3) +FPC*UB3 
MB (1) =QSD*CSL 

MB (2) =QSD*CM+FPC*UB3* (XBOD-XCG) 

MB(3)=QSD*CSN-FPC*OB2* (XBOD-XCG) 

GEES=-FB(1) / (MASS*G) 

IF(ABS(GEES) . GT . ABS (GMAX) )  GMAX=GEES 
C 

C . PRINT  TRAJECTORY  DATA 


00018500 

00018600 

00018700 

00018800 

00018900 

00019000 

00019100 

00019200 

00019300 

00019400 

00019500 

00019600 

00019700 

00019800 

00019900 

00020000 

00020100 

00020200 

00020300 

00020400 

00020500 

00020600 

00020700 

00020800 

00020900 

00021000 

00021100 

00021200 

00021300 

00021400 

00021500 

00021600 

00021700 

00021800 

00021900 

00022000 

00022100 

00022200 

00022300 

00022400 

00022500 

00022600 

00022700 

00022800 

00022900 

00023000 

00023100 

00023200 

00023300 

00023400 

00023500 

00023600 

00023700 

00023800 

00023900 

00024000 

00024100 

00024200 

00024300 

00024400 

00024500 

00024600 

00024700 

00024800 

00024900 

00025000 

00025100 

00025200 

00025300 

00025400 

00025500 

00025600 

00025700 
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c 

IF(T.LT.TTR)  GO  TO  14- 
TPR=TPR+DTPR 
30  CONTINUE 
H=-XE (3) 

VPE=SQRT (OE (1) **2+UE (2) **2+UE(3) **2) 

QDYN=0 . 5*RHO*VP*VP 
BXY=SQRT (B (1, 1) **2+B ( 1 , 2 ) **2) 

THETft=57 .295*ATAN2 (-B(l,3) ,BXY) 

ALPHAD=57 . 295*ATAN2 (SALP, CALP) 

UXY=SQRT (OE(l) **2+UE (2) **2) 

GAMMAD=57.295*ATAN2 {-UE(3) ,UXY) 

C  WRITE (12, 505)  T, H, XE ( 1) , XE (2 ) , VPE, VP, MACH, QDYN, GEES, 

C  $  GAMMAD, THETA, ALPHAD, CDS 

C  WRITE(6,203)  T, H, XE { 1 ) , XE (2 ) , VPE, VP, MACH, QDYN, GEES , 

C  $  GAMMAD, THETA, ALPHAD, CDS 

IF(IEND.NE.O)  GO  TO  31 
14  CONTINUE 
C 

C . EULER  ROTATION  FUNCTION  FOR  DIRECTION  COSINE  PROPAGATION 

C 

W2=WB (1) **2+WB (2) **2+WB (3) **2 
W=SQRT(W2) 

COSWT=COS (W*DT) 

SINWT=SIN(W*DT) 

COSWTM=l.-COSWT 
IF{W2 .GT. 1 .E-12)  GO  TO  22 
W2=1.E-12 
W=l.E-6 
22  CONTINUE 
C 

C . ANGULAR  MOMENTUM  CROSS  PRODUCT  TERMS 

C 

DO  20  K=l,3 

HB(K)=IN(K,  1)  *WB(1)+IN(K,2)*WB(2)+IN{K,  3)*WB{3) 

20  CONTINUE 

DO  21  1=1,3 
I1=E(I+1) 

I2=E(I+2) 

TEMP { I ) =WB (I 1 ) *HB ( I 2 ) -WB ( 12 ) *HB { I 1 ) 

21  CONTINUE 
C 

C . FORCE  RESOLUTIONS  TO  EULER  SYSTEM 

C  TRANSLATIONAL  ACCELERATIONS  AND  DERECTION  COSINE  ROTATION 
C 

D016  1=1,3 

FE(I)=FB(1)  *B(1,I)+FB(2)*B(2,I)+FB(3)*B(3,I) 

UEDOT(I)=FE(I) /MASS 
DO  17  J=l,3 
BN(I, J)=B(I, J) 

J1=E ( J+1) 

J2=E ( J+2) 

BDOT ( I , J) =DEL ( I , J) *COSWT+WB ( I ) *WB ( J) *COSWTM/W2+ 

S  (WB( Jl) *DEL(I, J2) -WB(J2) *DEL(I, Jl) ) *SINWT/W 
17  CONTINUE 
C 

C . ANGULAR  ACCELERATIONS  IN  BODY  AXES 

C 

WBDOT (I)=JN (I, 1) * (MB (1) -TEMP (1) ) 

$+JN{I,  2) * (MB (2) -TEMP (2) ) 

$+JN(I, 3) * (MB (3) -TEMP (3) ) 

16  CONTINUE 
C 

C . INTEGRALS 

C 

T=T+DT 
DO  11  1=1,3 

XE(I)=XE (I) +DT* (UE (I) +0.5*DT*OEDOT(I) ) 

UE ( I ) =UE ( I ) +DT*UEDOT ( I ) 

WB ( I ) =WB ( I ) +DT*WBDOT ( I ) 

DO  11  J=l,3 

B(I,  J)=BDOT(I,  1)  *BN(1,  J)+BDOT(I,2)*BN(2,  J)+BDOT(I,3)  *BN(3,  J) 
11  CONTINUE 


00025800 

00025900 

00026000 

00026100 

00026200 

00026300 

00026400 

00026500 

00026600 

00026700 

00026800 

00026900 

00027000 
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IF(T.GE.TEND)  GO  TO  15 
C  __ 

C . END  ACCELERATION  LOOP  -  TEST  FOR  GROUND  IMPACT 

C 

H=-XE(3) 

IF(H.GT.HMIN)  GO  TO  1 
15  CONTINUE 
IEND=1 
GO  TO  30 
31  CONTINUE 
T=-999. 

WRITE (12, 505)  T,GMAX 
505  FORMAT (13E10. 4) 

98  RETURN 
END 

SUBROUTINE  AERO 

COMMON/INERT/MASS,XCG,XBOD, IN (3, 3) , JN (3, 3) 

COMMON/BURN/  T 

COMMON/AEROS/  CBAR, SAREA, MACH 

COMMON/AEROV/  VPO, SALP, CALP, SBET, CBET, PB, QB,  RB 
COMMON/AEROF/  CN, CY, CA, CSL, CM, CSN, CDS 
COMMON/AEROCF/  CNA, CNA2, CNQ, CYB, CYB2 ,  CYR,  CAO,  CAA2 
COMMON/AEROCM/  CLO, CLP, CMO, CMA, CMA2, CMQ, CNB,  CNB2,  CNR 
COMMON/CHUTE/  IPTS, PCDS (50) , PT(50) 

COMMON/FBDRAG/  MPTS , PCDF ( 50 ) , PM(50) 

REAL  MASS, IN, JN, MACH 
C 

C . AERO  VARIABLES 

C 

SAC=SALP*CALP 
SAS=SALP*ABS (SALP) 

SBC=SBET*CBET 
SBS=SBET*ABS (SBET) 

RAD=CBAR/(2.*VP0) 

CNA2=0 . 

CNQ=0 . 

CYB2=0. 

CYR=0 . 

CM0=0 . 

CMA2=0 . 

CNB2=0 . 

C 

C . FOREBODY  ZERO-LIFT  DRAG  COEFFICIENT 

C 

1=0 

6  1=1+1 
IP=I+1 

IF(I.EQ.MPTS)  GO  TO  5 
IF(MACH.GT.PM(IP) )  GO  TO  6 

CA0=PCDF(I) + (PCDF(IP) -PCDF (I) ) * (MACH-PM(I) ) /(PM(IP) -PM (I) ) 
GO  TO  4 

5  CA0=PCDF{MPTS) 

4  CONTINUE 
C 

C . FORCE  AND  MOMENT  COEFFICIENTS 

C 

CN=CNA*SAC+CNA2*SAS+CNQ*QB*RAD 
CY=CYB*SBC+CYB2*SBS+CYR*RB*RAD 
CA=CA0+CAA2* ( 1 . -CALP**2*CBET**2) 

CSL=CL0+CLP*PB*RAD 

CM=CM0+CMA*SAC+CMA2*SAS+CMQ*QB*RAD 

CSN=CNB*SBC+CNB2*SBS+CNR*RB*RAD 

C 

C . PARACHUTE  DRAG-AREA 

C 

CDS=0. 

IF(T.LT.PT(1) )  GO  TO  1 
1  =  0 

3  1=1+1 
IP=I+1 

IF(I.EQ.IPTS)  GO  TO  2 
IF(T.GT.PT(IP) )  GO  TO  3 

CDS=PCDS (I)+(PCDS(IP) -PCDS(I) )*(T-PT(I) )/(PT(IP) -PT(I) ) 
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o  o  o 


WRITE  (*,*)•  T 


I 


IP 


CDS 


T<PT (IP) 


GO  TO  1 

2  CDS=PCDS (IPTS) 

1  CONTINUE 
RETURN 
END 

SUBROUTINE  DENS 

...DENSITY  AND  SPEED  OF  SOUND  VS.  ALTITUDE  -  DOMMASCH  EQUATION 
COMMON/REF/  H, SOUND, RHO, RHOZ 

RHO=RHOZ*EXP(-H/23111.+.294*SIN(H/28860. )+.213*SIN(H/86580. ) ) 

IF(H.GT.O)  GO  TO  7 

SOUND=1116.44 

RETURN 

7  IF(H.GT. 36152)  GO  TO  1 
T=518.688-(3.56616E-03)*H 
SOUND=49.02118*SQRT(T) 

RETURN 

1  IF(H.GT. 82345)  GO  TO  2 
SOUND=968 .08 

RETURN 

2  IF(H.GT. 155348)  GO  TO  3 
T=254 . 988+1 . 64592E-03*H 
SOUND=49.02118*SQRT(T) 

RETURN 

3  IF(H.GT. 175346)  GO  TO  4 
SOUND=1105.7 

RETURN 

4  IF(H.GT. 262448)  GO  TO  5 
T=988.088-2.46888E-03*H 
SOUND=49 . 02118+SQRT (T) 

RETURN 

5  IF(H.GT. 299516)  GO  TO  6 
SOUND=846.9 

RETURN 

6  T=-349. 812+2. 19456E-03+H 
SOUND=49 . 02118+SQRT (T) 

RETURN 

END 

FUNCTION  GRAV{H) 

RE=20855531.5 

GRAV=32 . 1741* (RE/ (H+RE) ) **2 

RETURN 

END 
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Input  Parameters  A 


★*★*★*★★*★ 

**********  SYSTEM  PARAMETERS 


88.18  SYSTEM  WEIGHT  (LBS) 

1.3125  FOREBODY  C.G.  (FT) 

2.625  FOREBODY  LENGTH  (FT) 

.1475  ROLL  INERTIA  IXX  (SLUG-FT**2) 

1.584  PITCH  INERTIA  lYY  (SLOG-FT**2) 

1.584  YAW  INERTIA  IZZ  (SLUG-FT**2) 


WEIGHT 

XCG 

XBOD 

IN(1, 1) 

IN(2,2) 

IN(3,3) 


**********  INITIAL  CONDITIONS 


328.1 
0.000 
484.65 
-85.46 
15.00 
0.000 
0.000 
0.000 

********* 


ALTITUDE  (FT)  ALT 

MINIMUM  ALTITUDE  OR  GROUND  LEVEL  (FT)  HMIN 

HORIZONTAL  VELOCITY  (FPS)  UE ( 1 ) 

EJECTION  VELOCITY  (FPS)  POSITIVE  DOWN  DE(3) 

PITCH  ANGLE  (DEG)  NOSE  UP  POSITIVE  THETA 

HEAD(+)  OR  TAIL(-)  WIND  (FPS)  VWIND(l) 

CROSSWIND  (FPS)  VWIND(2) 

DENSITY  (0.  FOR  STD.  ATMS.)  (SLUG/FT**3) DENSITY 


*  PROGRAMM  CONSTANTS 


**•***★★★★*★*★★★****★ 


0.0005 

0.10 

9.0 

******* 

******* 

******* 

******* 


INTEGRATION  TIME  STEP  (SEC) 
PRINT  INTERVAL  (SEC) 

MAXIMUM  TIME  OF  FLIGHT  (SEC) 


DT 

DTPR 

TEND 


***  FOREBODY  AERODYNAMIC  COEFFICIENTS  ******************** 


.6562 

REFERENCE  LENGTH  (FT)  Durchmesser 

CBAR 

.3382 

REFERENCE  AREA  (SQ.FT.)  Querschnitt 

SAREA 

2.780 

NORMAL  FORCE  CN-ALPHA  (/RAD) 

CNA 

0.000 

SIDE  FORCE  CY-BETA  (/RAD) 

CYB 

0.000 

AXIAL  FORCE  CA-ALPHA**2  (/RAD**2) 

CAA2 

0.000 

ROLL  TORQUE  COEFFICIENT  (NONDIM. ) 

CLO 

0.000 

ROLL  DAMPING  COEFFICIENT  (/RAD) 

CLP 

1.11 

PITCH  MOMENT  CM-ALPHA  (/RAD) 

CMA 

-10.0 

PITCH  DAMPING  (/RAD) 

CMQ 

0.000 

YAW  MOMENT  CN-BETA  (/RAD) 

CNB 

0.000 

YAW  DAMPING  (/RAD) 

CNR 

************* 

********** 

**********  fOREBODY  DRAG  VS  MACH  NUMBER 
********** 


NUMBER  OF  TABLE  INPUTS:  2 
NO.  MACH  NO.  DRAG  COEFFICIENT 

1  0.00  1.00 

2  1.00  1.00 

**********  parachute  DRAG-AREA  VS  TIME  ******************** 


DEPLOYMENT  TIME(SEC):  0.25 

NUMBER  OF  TABLE  INPUTS:  2  NOTE:  DRAG-AREA  VS.  TIME  IS 

NO.  TIME (SEC)  DRAG-AREA (SQ. FT. )  RELATIVE  TO  START  OF  DEPLOYMENT. 

1  0.00  0.3382 

2  0.10  10.333 
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Input  Parameters  B 


**********  SYSTEM  PARAMETERS  ******************** 


****★★★+★★ 

•*•*•****★*•**★★**★*★*★*★*'*•*•  +  *★*****★**★★★'( 
500.000 
4.000 
10.000 
5.000 
100.000 
100.000 


SYSTEM  WEIGHT  (LBS) 

FOREBODY  C.G.  (FT) 

FOREBODY  LENGTH  (FT) 

ROLL  INERTIA  IXX  (SLUG-FT**2) 
PITCH  INERTIA  lYY  (SLUG-FT**2) 
YAW  INERTIA  IZZ  (SmG-FT**2) 
************************************* 
****★*■**•** 

**********  initial  CONDITIONS 


WEIGHT 

XCG 

XBOD 

IN(1,1) 

IN(2,2) 

IN(3, 3) 


100.000 
0.000 
500.000 
0.000 
0.000 
0.000 
0.000 
0.000 

**********  proGRAMM  CONSTANTS  ******************** 

0.005  INTEGRATION  TIME  STEP  (SEC)  DT 

0.500  PRINT  INTERVAL  (SEC)  DTPR 

10.000  MAXIMUM  TIME  OF  FLIGHT  (SEC)  TEND 

**********  forebody  AERODYNAMIC  COEFFICIENTS  ******************** 


ALTITUDE  (FT)  ALT 

MINIMUM  ALTITUDE  OR  GROUND  LEVEL  (FT)  HMIN 

HORIZONTAL  VELOCITY  (FPS)  UE ( 1 ) 

EJECTION  VELOCITY  (FPS)  POSITIVE  DOWN  UE(3) 

PITCH  ANGLE  (DEG)  NOSE  UP  POSITIVE  THETA 

HEAD(+)  OR  TAIL(-)  WIND  (FPS)  VWIND(l) 

CROSSWIND  (FPS)  VWIND(2) 

DENSITY  (0.  FOR  STD.  ATMS.)  (SLUG/FT**3) DENSITY 


10.000 

REFERENCE  LENGTH  (FT) 

CBAR 

1.000 

REFERENCE  AREA  (SQ.FT.) 

SAREA 

0.000 

NORMAL  FORCE  CN-ALPHA  (/RAD) 

CNA 

0.000 

SIDE  FORCE  CY-BETA  (/RAD) 

CYB 

0.000 

AXIAL  FORCE  CA-ALPHA**2  (/RAD*»2) 

CAA2 

0.000 

ROLL  TORQUE  COEFFICIENT  (NONDIM. ) 

CLO 

0.000 

ROLL  DAMPING  COEFFICIENT  (/RAD) 

CLP 

-2.000 

PITCH  MOMENT  CM-ALPHA  (/RAD) 

CMA 

-200.000 

PITCH  DAMPING  (/RAD) 

CMQ 

0.000 

YAW  MOMENT  CN-BETA  (/RAD) 

CNB 

0.000 

YAW  DAMPING  (/RAD) 

CNR 

******************************************************** 

★  +  ★★***★★*★****★***■*• 

**********  fOREBODY  drag  VS  MACH  NUMBER  ******************** 

**********  ******************** 

********************************************************************** 
NUMBER  OF  TABLE  INPUTS:  2 
NO.  MACH  NO.  DRAG  COEFFICIENT 

1  0.00  0.09 

2  2.00  0.09 

********************************************************************** 
**********  ******************** 

**********  parachute  drag-area  vs  time  ******************** 

**********  ******************** 

********************************************************************** 
DEPLOYMENT  TIME(SEC):  0.20 
NUMBER  OF  TABLE  INPUTS:  2 
NO.  TIME (SEC)  DRAG-AREA(SQ.FT. ) 

1  0.00  0.20 

2  1.00  500.00 


NOTE:  DRAG-AREA  VS.  TIME  IS 
RELATIVE  TO  START  OF  DEPLOYMENT. 
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APPENDIX  B 

"MODSIM  III  Translation  of  Purvis  6-DOF  Model 


MAIN  MODULE  sixDOF; 

FROM  globalMod  IMPORT  jumper; 

FROM  globalMod  IMPORT  initializeData; 

BEGIN 


initializeData; 

NEW  (jumper) ; 

ASK  jumper  TO  jump; 

END  (MAIN)  MODULE  {6D0F}. 
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DEFINITION  MODULE  globalMod; 

FROM  jumperMod  IMPORT  jumperObj; 
CONST 


re  =  20855531.5; 


TYPE 


VAR 


eType  = 
delType  = 
matrixType  = 
vectorType  = 


ARRAY  INTEGER  OF  INTEGER; 

ARRAY  INTEGER,  INTEGER  OF  REAL; 
ARRAY  INTEGER,  INTEGER  OF  REAL; 
ARRAY  INTEGER  OF  REAL; 


i,  j  :  INTEGER; 

jumper  :  jumperObj; 

e  :  eType ; 

del  :  delType; 


PROCEDURE  initializeData; 


END  {DEFINITION}  MODULE  {globalMod}. 
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IMPLEMENTATION  MODULE  globalMod; 

PROCEDURE  initializeData; 

BEGIN 

NEW  (e,  1..5); 

NEW  (del,  1..3, 

e[l]  :=  1; 
e[2]  :=  2; 
e[3]  :=  3; 
e[4]  :=  1; 
e[5]  :=  2; 

FOR  i  :=  1  TO  3 

FOR  j  :=  1  TO  3 
IF  i  =  j 

del[i,j]  :=  1.0; 

ELSE 

del[i,j]  :=  0.0; 

END  IF; 

END  FOR; 

END  FOR; 

END  PROCEDURE  (initializeData); 

END  (IMPLEMENTATION)  MODULE  (globalMod). 
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DEFINITION  MODULE  calcMod; 

PROCEDURE  gravCalc  TiN  a  :  REAL)  :  REAL; 

PROCEDURE  densityCalc  (IN  h,  rhoz  ;  REAL;  OUT  rho,  sound  :  REAL) 
END  (DEFINITION)  MODULE  (calcMod). 


IMPLEMENTATION  MODULE  calcMod; 


FROM  MathMod  IMPORT  POWER,  SIN,  COS,  SQRT,  EXP; 

FROM  globalMod  IMPORT  re; 

PROCEDURE  gravCalc  (IN  a:  REAL)  :  REAL; 

BEGIN 

RETURN  32.1741*POWER(re/(a+re)  ,  2.0); 

END  PROCEDURE  (gravCalc); 

PROCEDURE  densityCalc  (IN  h,  rhoz  :  REAL;  OUT  rlio,  sound  :  REAL)  ; 
VAR 


t  :  REAL; 


BEGIN 

rlio  :=  rhoz  *  EXP  (-1 . 0*h/23111 .0  +  0.294  *  SIN  (h/288  60 . 0 )  +  0.213  * 
SIN(h/86580.0) ) ; 

IF  h  >  0.0 

t  :=  518.688  -  (3 . 56616E-03) ‘h; 
sound  :=  49.02118  *  SQRT(t); 

IF  h  >  36152.0 

sound  :=  968.08; 

IF  h  >  82345.0 

t  :=  254.988  +  (1 . 64592E-03) ♦h; 
sound  :=  49.02118  *  SQRT(t); 

IF  h  >  155348.0 

sound  :=  1105.0; 

IF  h  >  262448.0 

sound  :=  846.9; 

IF  h  >  299516.0 

t  :=  -349.812  +  2 . 19456E-03*h 
END  IF; 

END  IF; 

END  IF; 

END  IF; 

END  IF; 

ELSE 

sound  ;=  1116.44; 

END  IF; 

END  PROCEDURE  (densityCalc); 

END  (IMPLEMENTATION)  MODULE  (calcMod) . 
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DEFINITION  MODULE  jumperMod; 

FROM  globalMod  IMPORT  matrixType,  vectorType; 
TYPE 

jumperObj  =  OBJECT 

11, 

12, 
lend, 
ipts, 
ip» 
jl- 
j2, 
k, 

loop, 

mpts, 

test  :  INTEGER; 

alt, 

alpha, 

alphad, 

beta, 

bxy, 

ca, 

cao, 

caa2, 

calpha, 

cbar, 

cbeta, 

cby2 , 

cds, 

do, 

dp, 

cm, 
cma, 
cma2 , 
cmo, 
cmq, 

cn, 

cna, 
cna2 , 

cnb, 
cnb2, 

cnq, 

cnr, 
coswt, 
coswtm, 

CP/ 

cphi, 

csl, 

csn, 

ct, 

cy, 

cyb, 

cyb2, 

cyr, 

dens, 

deptime, 

dt, 

dtpr, 

fpc, 

g/ 

qaminad , 
gees, 
gmax, 
h, 

hmin, 

mach, 

mass, 

pb, 

phi. 


psi, 

q- 

qb, 

qdyn, 

qs, 

qsd, 
rad, 
rb, 

rho, 
rhoz, 
sac, 
sas, 
sarea, 
salpha, 
sbc, 
sbeta, 
sbs, 
sinwt, 
sound, 
sp, 
phi, 

St., 

t, 

tend, 

theta, 

tar, 

ubl, 

ub2, 

ub3, 

uel, 

ue2, 

ue3, 

uxy, 

vp, 

vpl3, 

vpe, 

vpo, 

w, 

w2, 

weight, 

xbod, 

xcg  :  REAL; 

pcds, 

Pt, 

pcdf, 

pm  :  vectorType;  {1X2} 

xe, 

ue, 

wb, 

vwind, 

temp, 

fb, 

m, 

mb, 

fe, 

uedot, 

wbdot, 

hb  :  vectorType;  {1X3} 

in, 

jn, 

b, 

bn, 

bdot  :  matrixType;  {3X3} 
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ASK  METHOD  Objinit;  END  OBJECT  { jumperObj ) ; 

ASK  METHOD_  jump; 

END  {DEFINITION}  MODULE  { jumperMod} . 
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IMPLEMENTATION  MODOLE  jumperMod; 

FROM  MathMod  IMPORT  EXP,  SIN,  COS,  POWER,  SQRT,  ATAN2; 
FROM  globalMod  IMPORT  re,  e,  del,  i,  j; 

FROM  calcMod  IMPORT  gravCalc,  densityCalc; 

OBJECT  jumperObj; 

ASK  METHOD  Objinit; 

BEGIN 


NEW(pcdf,  1..2); 
NEW (pm  ,  1..2); 

NEW ( f b  ,  1 . . 3 ) 
NEW{fe  ,  1..3) 
NEW(hb  ,  1..3) 
NEW (mb  ,  1..3) 
NEW (temp  ,  1..3) 
NEW(ue  ,  1..3) 
NEW(uedot,  1..3) 
NEW (vwind,  1 . . 3) 
NEW(wb  ,  1..3) 
NEW(wbdot,  1..3) 
NEW(xe  ,  1..3) 


NEW(pcds,  1. 

.4)  ; 

NEW(pt  ,  1. 

.4)  ; 

NEW (in  ,  1. 

.3,  1..3); 

NEW{jn  ,  1. 

.3,  1..3); 

NEW(b  ,  1. 

.3,  1..3); 

NEW (bn  ,  1. 

.3,  1..3); 

NEW(bdot,  1. 

.3,  1..3); 

(  system  inertial  properties  } 


weight  : 

=  88.18; 

{ 

parachute-payload  system  weight 

mass  : 

=  weight/32. 

17; 

xcg  : 

=  1.3125; 

{ 

forebody  c.g.  (ft)  ) 

xbod  : 

=  2.625; 

{ 

forebody  length  (ft)  ) 

FOR  i  := 

1  TO  3 

FOR 

j  :=  1  TO  3 

jn{i,j]  := 

0.0; 

in(i,j]  := 

0.0; 

END 

(j)  for; 

END  (i) 

FOR; 

in(l,l] 

:=  0.1475;  { 

roll 

inertia  Ixx  (slug-ft''2) 

) 

in(2,2] 

:=  1.584;  { 

pitch 

inertia  lyy  {slug-ft''2) 

) 

in[3,3] 

:=  1.584;  { 

yaw 

inertia  Izz  {slug-ft''2) 

} 

FOR  i  :=  1  TO  3 

jn[i,i]  :=  1.0  /  in  [i,i]; 
END  FOR; 


{  initial  conditions  ) 


FOR  i  1  TO  3 


xe[i]  := 

0.0; 

{ 

(ft) 

1 

ue[il  := 

0.0; 

{ 

(fps) 

1 

ejection  velocity  positive 

down) 

wb [ i ]  :  = 

0.0; 

{ 

???  } 

vwind [i]  := 

0.0; 

{ 

(fps) 

1 

????  ) 

END  FOR; 


down  range,  2:  off  range,  3:  altitude  loss) 
horizontal  velocity,  2:  lateral  velocity,  3 

head  {+)  or  tail  (-)  wind,  2:  crosswind,  3: 


alt 
hmin 
ue[l] 
ue  [3] 
theta 


=  328.1; 

=  0.0; 

=  484.65; 
=  -85.46; 
=  15.00; 


{altitude  (ft)  ) 

{ground  level  (ft)  > 

{pitch  angle  (deg)  nose  up  positive) 
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vwind I 1 ]  : =  0.0; 
vwind[2]  „:=  0.0; 

dens  :=  0.0;  {density  (0  for  standard  atms)  in  slug/ft^3} 

rhoz  :=  0.002378;  (????) 

IF  dens  <>  0.0 

rhoz  :=  dens  *  EXP  (alt/23111.0  -  0.295  *  SIN (alt/28860 . 0)  -  0.213  * 
SIN{alt/86580.0) ) 

END  IF; 

{  program  constants  } 

dt  :=  0.0005;  {  integration  time  step  (sec)  } 

dtpr  :=  0.1;  {  print  interval  (sec)  ) 

tend  :=  9.0;  {  max  time  for  flight  (sec)  } 

{  forebody  aerodynamic  coefficients  } 

cbar  :=  0.6562;  {  reference  length  (ft)  } 

sarea  :=  0.3382;  {  reference  area  (ft''2)  ) 

cna  :=  2.78;  {  normal  force  cn-alpha  (/rad)  ) 

cyb  :=  0.0;  {  side  force  cy-beta  (/rad)  } 

caa2  :=  0.0;  {  axial  force  ca-alpha''2  (/rad^2)  ) 

do  :=  0.0;  {  roll  torque  coefficient  (dimensionless)  } 

clp  :=  0.0;  {  roll  damping  coefficient  (/rad)  ) 

cma  :=  1.11;  {  pitch  moment  cm-alpha  (/rad)  } 

cmq  :=  -10.0;  {  pitch  damping  (/rad)  } 

cnb  :=  0.0;  {  yaw  moment  cn-beta  (/rad)  } 

cnr  :=  0.0;  {  yaw  damping  (/rad)  ) 

{  forebody  drag  versus  mach  number  table  ) 

mpts  :=  2; 

pm[l]  :=  0.00;  (  mach  number  } 

pm[2]  :=  1.00; 

pcdfd]  :=  1.00;  (  drag  coefficient  ) 

pcdf[21  :=  1.00; 

{  parachute  drag-area  versus  time  table  ) 

deptime  :=  0.25;  {  deployment  time  ) 
ipts  ;=  2; 

IF  deptime  <  0.0 

deptime  :=  0.0; 

ipts  :=  ipts  +2;  {  ipts  =  4  } 

pt[l]  ;=  0.0; 

pt[2]  :=  deptime; 

pt[31  :=  0.00  +  deptime; 

pt(4]  :=  0.10  +  deptime; 

pcds[l]  :=  0.0; 

pcds[2]  :=  0.0; 

pcds(3]  :=  0.3382; 

pcds[4]  :=  10.333; 

IF  deptime  <=  0.0 

pt[l]  :=  0.0  +  deptime; 
pt[2]  :=  0.10  +  deptime; 
pcds[l]  :=  0.3382; 
pcds[2]  :=  10.333; 

END  IF; 

END  IF; 

IF  deptime  >  0.0 

ipts  :=  ipts  +2;  {  ipts  =  4  } 

pt(l]  :=  0.0; 

pt[2]  :=  deptime; 

pt[3]  :=  0.00  +  deptime; 

pt[4]  :=  0.10  +  deptime; 

pcds [ 1 ]  : =  0.0; 

pcds[2]  :=  0.0; 

pcds[3]  :=  0.3382; 

pcds [4]  :=  10.333; 

END  IF; 
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(  convert  EULER  ANGLES  to  di-rection  cosines  ) 


psi 

- 

0.0; 

phi 

= 

0.0; 

rad 

1.0/57.295; 

St. 

=  SIN 

{theta*rad) ; 

ct 

- 

COS 

(theta*rad) ; 

sp 

= 

SIN 

(psi*rad) ; 

cp 

COS 

(psi*rad) ; 

phi  : 

= 

SIN  (phi*rad); 

cphi 

= 

COS 

(phi* rad) ; 

xe  [3] 

= 

-1.0 

*  alt; 

b[l,l] 

= 

cp  * 

ct; 

b[l,2] 

sp  * 

ct; 

b[l,3] 

= 

-1.0 

*  St.; 

b[2,l] 

-1.0 

*  sp  *  cphi  +  cp  *  St . 

*  phi 

b[2,2] 

= 

cp  * 

cphi  +  sp  *  St.  *  phi; 

b[2,31 

= 

ct  * 

phi; 

b[3,  1] 

= 

sp  * 

phi  +  cp  *  St.  *  cphi; 

b[3,2] 

= 

-1.0 

*  cp  *  phi  +  sp  *  St. 

*  cphi 

b[3,3] 

= 

ct  * 

cphi ; 

END  METHOD  {Objinit}; 


ASK  METHOD  jump; 
BEGIN 


gees 

cds 

t 

gmax 

lend 

tar 

tend 

alt 

h 

test 


0.0; 

0.0; 

0.0; 

0.0; 

0; 

t; 

tend  -  0.1  *  dt; 

-1.0  *  xe[3];  {  alt  =  328.1  ft  ) 
alt; 

0; 


OUTPUT  (tend,  "  ",  hmin,  "  ",  alt); 

OUTPUT  (ue(l],  "  ",  ue(3],  "  ",  vwindt2]); 


(BEGIN  TRAJECTORY  LOOP} 


WHILE  test  =  0 


IF  h  >  hmin 


g  :=  gravCalc(h); 

densityCalc  (h,  rhoz,  rho,  sound) ; 

pb  : =  wb ( 1 ] ; 

qb  :=  wb[2]; 

rb  :  =  wb  ( 3  ]  ; 

uel  :=  ue[l]  -  vwind[l]; 

ue2  :=  ue[2]  -  vwind[2]; 

ue3  :=  ue[3]  -  vwind(3); 

vp  :=  SQRT(POWER(uel,2.0)  +  POWER(ue2, 2 . 0)  " 

vpo  :=  vp; 

mach  :=  vp /sound ; 

ubl  ;=  b[l,l]*uel  +  b[l,2]*ue2  +  b[l,3]*ue3; 
ub2  :=  b[2,l]*uel  +  b[2,2]*ue2  +  b[2,3]*ue3; 
ub3  :=  b[3,l]*uel  +  b[3,2]*ue2  +  b[3,3]*ue3; 
vpl3  :=  SORT (POWER(ubl, 2.0)  +  POWER (ub3 , 2 . 0 ) ) ; 

(  USE  SIN (ALPHA)  for  ALPHA  and  COS (BETA)  for  BETA  ) 

IF  vpo  <  l.OE-06 

vpo  :=  l.OE-06; 

END  IF; 

sbeta  :=  ub2  /  vpo; 


POWER(ue3,2.0) ) ; 
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cbeta  :=  vpl3  /  vpo; 
beta  :=  sbeta; 


IF  vpl3  <  l.OE-06 

vpl3  ;=  l.OE-06; 

END  IF; 

salpha  :=  ub3  /  vpl3; 
calpha  :=  ubl  /  vpl3; 
alpha  :=  salpha; 

{  AERODYNAMIC  and  BODY  FORCES  AND  MOMENTS  } 

(  ISOLATED  BODY  AERODYNAMICS  } 

{  AERO  ROUTINE  ) 

sac  :=  salpha  *  calpha; 

sas  :=  salpha  *  ABS  (salpha); 

sbc  :=  sbeta  *  cbeta; 

sbs  :=  sbeta  *  ABS (sbeta) ; 

rad  :=  cbar  /  (2.0*vpo); 

cna2  :=  0.0; 

cnq  : =  0.0; 

cyb2  : =  0.0; 

cyr  : =  0.0; 

cmo  : =  0.0; 

cma2  : =  0.0; 

cnb2  ; =  0 . 0 ; 

(  FOREBODY  XERO-LIFT  DRAG  COEFFIECIENT  ) 

i  :=  0; 
loop  :=  0; 

WHILE  loop  =  0 

i  :=  i  +  l.- 

ip  :=  i  +  1; 

IF  i  =  mpts 

cao  :=  pcdf(2];  {  when  i  =  ipts  ) 
loop  :=  1; 

ELSE 

IF  mach  <=  pm(ip) 

cao  :=  pcdf [i] + (pcdf [ip] -pcdf [i) ) * (mach- 

pm[i] ) / (pm[ip] -pm[i] ) ; 

loop  :=  1; 

END  IF; 

END  IF; 

END  WHILE; 

cn  :=  cna  *  sac  +  cna2  *  sas  +  cnq  *  qb  *  rad; 

cy  :=  cyb  +  sbc  +  cyb2  *  sbs  +  cyr  *  rb  *  rad; 

ca  :=  cao  +  caa2  *  (1.0  -  POWER (calpha, 2 . 0)  *  POWER (cbeta, 2 . 0) ) 
csl  :=  do  +  clp  *  pb  *  rad; 

cm  :=  cmo  +  cma  *  sac  +  cma2  *  sas  +  cmq  *  qb  *  rad; 

csn  :=  cnb  *  sbc  +  cnb2  *  sbs  +  cnr  *  rb  *  rad; 

{  PARACHUTE  DRAG-AREA  ) 

cds  : =  0.0; 
i  :=  0; 
loop  :=  0; 

WHILE  loop  =0 

IF  t  <  pt[l] 

loop  :=  1; 

ELSE 

i  ;=  i  +  1; 
ip  :=  i  +  1; 

IF  i  =  ipts 

cds  ;=  pcds[ipts]; 
loop  :=  1; 

ELSE 
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IF  t  <=  pt  [ip] 

_  ■  cds  :=  pcds[i]  +  {pcds [ip] -pcds [i] )  *  {t-pt[i])  / 

(pt [ip] -pt  [i] )  ; 

loop  :=  1; 

END  IF; 

END  IF; 

END  IF; 

END  WHILE; 

{  END  AERO  ROUTINE  ) 

q  :=  0.5  *  rho  *  vp; 
fpc  :=  -I.O  *  q  *  cds; 

qs  :=  0.5  *  rho  *  POWER  (vp,  2.0)  *  sarea; 
qsd  :=  qs  *  char; 

fb[l]  ;=  -1.0  *  qs  *  ca  +  mass  *  g  *  b[l,3]  +  fpc  *  ubl; 

fb[2]  :=  qs  *  cy  +  mass  *  g  *  b[2,3]  +  fpc  *  ub2; 

fb[3]  :=  -1.0  *  qs  *  cn  +  mass  *  g  *  b[3,3)  +  fpc  *  ub3; 

inb[l]  :=  qsd  *  csl; 

mb[2]  :=  qsd  *  cm  +  fpc  ♦  ub3  *  (xbod  -  xcg) ; 

mb[3]  ;=  qsd  *  csn  -  fpc  ♦  ub2  ♦  (xbod  -  xcg) ; 
gees  :=  -1.0  *  fb[l]  /  (mass  *  g) ; 

IF  ABS  (gees)  >  ABS  (graax) 
gmax  :=  gees; 

END  IF; 

IF  lend  =  0 

IF  t  >=  tar 

tar  :=  tar  +  dtpr; 
h  : - 1.0  *  xe[3]; 

vpe  :=  SORT  (  POWER  (uel,2.0)  +  POWER  (ue2,2.0)  +  POWER 

(ue3,2.0)  ); 

qdyn  :=  0.5  *  rho  *  POWER  (vp,2.0); 

bxy  :=  SORT  (  POWER  (b[l,l),  2.0)  +  POWER  (b[l,2],2.0)  >; 
theta  :=  57.295  ♦  ATAN2  (  (-1.0  *  b[l,3]),  bxy  ); 
alphad  :=  57.295  *  ATAN2  (  salpha,  calpha  ); 
uxy  :=  SORT  (  POWER  (ue[l),2.0)  +  POWER  (ue[2],2.0)  ); 
gammad  :=  57.295  *  ATAN2  (  (-1.0  *  ue[3)),  uxy  ); 

OUTPUT  (t,  '■  ",  h,  "  ",  xe[l],  "  ",  xe[2],  "  ",  vpe,  "  ", 
vp,  "  ",  mach,  "  ",  qdyn,  "  ",  gees,  "  ",  gammad,  "  ",  theta,  "  ",  alphad,  "  ",  cds); 

END  IF; 

{  EULER  ROTATION  FUNCTION  FOR  DIRECTION  COSINE  PROPOAGATION  ) 

w2  :=  POWER(wb[l] ,2.0)  +  POWER(wb [2] , 2 . 0)  +  POWER (wb [3] , 2 . 0 ) ; 

w  :=  SORT  (w2); 

coswt  :=  COS  (w*dt); 

sinwt  :=  SIN  (w*dt) ; 

coswtm  ;=  I.O  -  coswt; 

IF  w2  <  l.OE-12 

w2  :=  l.OE-12; 
w  :=  l.OE-06; 

END  IF; 

[  ANGULAR  MOMENTUM  CROSS  PRODUCT  TERMS  } 

FOR  Ic  :=  1  TO  3 

hb  [k]  :=  in[k,l]  *  wb[l]  +  in[k,2]  *  wb[2]  +  in[k,3]  * 

wb[31  ; 

END  FOR; 

FOR  i  :=  1  TO  3 

11  :=  e[i+l]; 

12  :=  e[i+2]; 

temp[i]  ;=  wb[il]  *  hb[i2]  -  wb[i2]  *  hb[il]; 

END  FOR; 

(  FORCE  RESOLUTION  TO  EULER  SYSTEM  } 

(  TRANSLATIONAL  ACCELERATION  AND  DIRECTION  COSINE  ROTATION  ) 

FOR  i  :=  1  TO  3 
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b[2,i]  +  fb[3]  *  b[3,i]; 


fe[i]  :=  fb[l]  *  b[l,i]  +  fb[2]  * 

_  uedot[i]  :=  fe[il  /  mass; 

FOR  j  :=  1  TO  3 

bn[i, j]  :=  b[i, j]; 
jl  :=  e[j  +  l]; 
j2  :=  etj+2]; 

bdot[i,j]  :=  del [i, j ] *coswt  +  wb [i] *wb [ j ] *coswtm/w2  +  ( 
wb [ jl] *del [i, j2]  -  wb [ j2 ] *del [i, j 1]  )*  sinwt/w; 

END  FOR; 


{  ANGULAR  ACCELERATION  IN  BODY  AXES  } 


wbdot[i]  :=  jn [i, 1] * (mb [1] -temp [ 1] )  +  jn [i, 2 ] * (mb [2 ] - 
temp[2])  +  jn [i, 3] * (mb [3] -temp [3] ) ; 

END  FOR; 


bdot [i, 3] *bn [3, j ] ; 


(  INTEGRALS  } 
t  :=  t  +  dt; 

FOR  i  :=  1  TO  3 

xe(i]  :=  xe[i]  +  dt* (ue [i] +0 . 5*dt*uedot [i] ) ; 
ue[i]  :=  ue[i]  +  dt*uedot[i]; 
wb[i]  :=  wb(i]  +  dt*wbdot[i]; 

FOR  j  :=  1  TO  3 

b[i,j]  :=  bdot[i, 1] *bn[l, j ]  +  bdot [i, 2 ] *bn [2, j ]  + 

END  FOR; 

END  FOR; 


IF  t  >=  tend 

lend  :=  1; 

ELSE 

h  ;=  -1.0  *  xe(31; 

END  IF; 

ELSE  (  when  lend  =  1  } 
h  :=  -1.0  *  xe(3); 

vpe  :=  SORT  (  POWER  (uel,2.0)  +  POWER  (ue2,2.0)  +  POWER  (ue3,2.0) 

); 

qdyn  :=  0.5  *  rho  *  POWER  (vp,2.0); 

bxy  :=  SORT  (  POWER  (b(l,l],2.0)  +  POWER  (b[l,2],2.0)  ); 
theta  :=  57.295  ♦  ATAN2  (  (-1.0  *  b[l,3]),  bxy  ); 
alphad  :=  57.295  *  ATAN2  (  salpha,  calpha  ); 
uxy  :=  SORT  (  POWER(ue[lJ ,2.0)  +  POWER (ue [2] , 2 . 0)  ); 
gammad  :=  57.295  *  ATAN2  (  (-1.0  *  ue[3]),  uxy  ); 

OUTPUT  (t,  ”  ",  h,  "  ",  xe[l],  "  ",  xe[2],  "  ",  vpe,  "  ",  vp,  " 

",  mach,  "  ",  qdyn,  "  ",  gees,  "  ",  gammad,  "  ",  theta,  "  ",  alphad,  "  ",  cds); 

test  :=  1; 

END  IF  (lend); 

ELSE  (  when  h  <=  hmin  } 
test  :=  1; 

END  IF  (hmin); 


END  WHILE  (test); 

END  METHOD  (jump); 

END  OBJECT  (jumperObj); 

END  (IMPLEMENTATION)  MODULE  { jumperMod) . 
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APPENDIX  C 

Comparisons  of  Purvis  6-DOF  and 
MODSIM  6-DOF  State  Variable  Outputs 


Figure  20.  Altitude  Vs.  Time 
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Figure  22.  Altitude  Vs.  Down  Range 


Figure  23.  Descent  Velocity  Vs.  Time 
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MACH  NUMBER 


Figure  24.  Airspeed  Vs.  Time 
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Figure  26.  Dynamic  Pressure  Vs.  Time 


2  4  6  8  To 

TIME 

Figure  27.  Axial  Acceleration  Vs.  Time 
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APPENDIX  D 

MODSIM  III  Paratrooper  Object  Modules 

MAIN  MODULE  vortex; 

FROM  inputMod  IMPORT  readData,  disposeStreams; 

FROM  globalMod  IMPORT  i,  Numberof Planes,  nu,  knotconv,  initializeData,  repeat 

FROM  AirplaneMod  IMPORT  ClVObj; 

FROM  VortexMod  IMPORT  RightVortexOb j ,  Lef tVortexObj ; 

FROM  MathMod  IMPORT  pi; 

FROM  VortexControlMod  IMPORT  Airdrop; 

FROM  SimMod  IMPORT  ResetSimTime,  StartSimulation; 

FROM  UtilMod  IMPORT  DateTime; 

VAR 

BEGIN 

{  -  Start  the  input  questions  and  set  up  the  random  seeds  -  } 

readData; 

initializeData; 

FOR  repeat  :=  1  TO  50; 

ResetSimTime (0.0); 

{  -  Create  the  Vortex  Control  Object  named  Airdrop  -  ) 

NEW  (Airdrop) ; 

{  -  Schedule  the  first  event  to  intiate  the  simulation  -  ) 

TELL  Airdrop  TO  Fly; 

StartSimulation; 

DISPOSE  (Airdrop) ; 

END  FOR; 
disposeStreams; 

END  (MAIN)  MODULE  (Vortex) . 
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DEFINITION  MODULE  globalMod; 

FROM  RandMod  IMPORT  RandomObj; 

FROM  VortexMod  IMPORT  RightVortexOb j ,  LeftVortexObj ; 

FROM  right JumperMod  IMPORT  right JumperObj ; 

FROM  left JumperMod  IMPORT  left JumperObj ; 

CONST 


re  =  20855531.5; 
nu  =  0.0001654; 

knotconv  =  1.69085;  (Converts  knots  to  ft/sec} 


TYPE 


eType  =  ARRAY  INTEGER  OF  INTEGER; 

delType  =  ARRAY  INTEGER,  INTEGER  OF  REAL; 
matrixType  =  ARRAY  INTEGER,  INTEGER  OF  REAL; 
vectorType  =  ARRAY  INTEGER  OF  REAL; 

encounterType  =  RECORD 

airplane  ;  INTEGER; 
side  :  STRING; 

position  :  INTEGER; 

END  RECORD  (encounterType); 

ElementPositionType  =  RECORD 
ElementPosNum  :  INTEGER; 

Intrail  :  REAL; 

CrossTrack  ;  REAL; 

END  RECORD; 


ElementGeometryType 

Format ionPositionType 
PositonNumber 
Intrail 
CrossTrack 
END  RECORD; 


=  ARRAY  INTEGER  OF 

=  RECORD; 

INTEGER; 

REAL; 

REAL; 


ElementPositionType; 


FormationGeometryType  =  ARRAY  INTEGER  OF  FormationPositionType 


VAR 


NumberofPlanes  : 
PlanesPerElement  : 
NumberOfElements  : 
i,  j,  repeat  : 

Elementspacing  ; 
ElementGeometry  ; 
FormationGeometry  ; 

CrossWindl  : 
CrossWind2  : 
CrossWindS  : 
ShearAltl  : 
ShearAlt2  : 
StandDevl  : 
StandDev2  : 
StandDev3 

trailBox  : 
lateralBox  : 
altitudeBox  : 
Headwind  : 
rho  : 
altitude  : 
vfk  : 
weight  : 
RunLength  : 
vsl  : 
vs2  : 
vs3  : 


INTEGER; 

INTEGER; 

INTEGER; 

INTEGER; 

REAL; 

ElementGeometryType; 

FormationGeometryType; 

REAL; 

REAL; 

REAL; 

REAL; 

REAL; 

REAL; 

REAL; 

REAL; 

INTEGER; 

INTEGER; 

INTEGER; 

REAL; 

REAL; 

REAL; 

REAL; 

REAL; 

REAL; 

REAL; 

REAL; 

REAL; 
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rightJumper  ___ 

right JumperObj ; 

left Jumper 

left JumperObj  ; 

e 

eType; 

del 

delType; 

seedl 

RandomObj ; 

seed2 

RandomOb j ; 

seed3 

RandomObj ; 

seed4 

RandomObj ; 

windseedl 

RandomObj ; 

windseed2 

RandomObj ; 

windseed3 

RandomObj ; 

trailseed 

RandomObj ; 

lateralseed 

RandomObj ; 

timeseed 

RandomObj ; 

PROCEDURE  initializeData; 

END  (DEFINITION}  MODULE  {globalMod}. 
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IMPLEMENTATION  MODULE  globalMod; 


FROM  RandMod  IMPORT  FetchSeed; 

FROM  inputMod  IMPORT  jumperseed; 

PROCEDURE  initializeData; 

BEGIN 

NEW  (e,  1..5); 

NEW  (del,  1..3,  1..3); 

NEW  (seedl); 

NEW  (seed2); 

NEW  (seed3); 

NEW  (seed4); 

NEW  (windseedl ) ; 

NEW  (wind3eed2 ) ; 

NEW  (windseed3 ) ; 

NEW  (trailseed) ; 

NEW  (lateralseed) ; 

NEW  (timeseed) ; 

ASK  seedl  TO  SetSeed  (FetchSeed  (  jumperseed. Uniformint  (1,10)  ) ) ; 

ASK  seed2  TO  SetSeed  (FetchSeed  (  jumperseed. Uniformint  (1,10)  )); 

ASK  seed3  TO  SetSeed  (FetchSeed  (  jumperseed . Uniformint  (1,10)  ) ) ; 

ASK  seed4  TO  SetSeed  (FetchSeed  (  jumperseed . Uniformint  (1,10)  )  )  ; 

ASK  windseedl  TO  SetSeed  (FetchSeed  (  jumperseed. Uniformint  (1,10)  ) ) ; 

ASK  windseed2  TO  SetSeed  (FetchSeed  (  jumperseed. Uniformint  (1,10)  ) ) ; 

ASK  windseed3  TO  SetSeed  (FetchSeed  (  jumperseed. Uniformint  (1,10)  ) ) ; 

ASK  trailseed  TO  SetSeed  (FetchSeed  (  jumperseed . Uniformint  (1,10)  )); 

ASK  lateralseed  TO  SetSeed  (FetchSeed  (  jumperseed. Uniformint  (1,10)  )) 
ASK  timeseed  TO  SetSeed  (FetchSeed  (  jumperseed. Uniformint  (1,10)  )); 


e[l]  :=  1; 
e[2]  :=  2; 
e[3]  :=  3; 
e[4!  :=  1; 
e[5]  :=  2; 

FOR  i  :=  1  TO  3 

FOR  j  :=  1  TO  3 
IF  i  =  j 

del(i,j]  :=  1.0; 

ELSE 

del[i,j]  :=  0.0; 
END  IF; 

END  FOR; 

END  FOR; 

END  PROCEDURE  {initializeData}; 

END  {IMPLEMENTATION}  MODULE  {globalMod}. 
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DEFINITION  MODULE  calcMod; 


PROCEDURE  gravCalc  (IN  a  :  REAL)  :  REAL; 

PROCEDURE  densityCalc  (IN  h,  rhoz  :  REAL;  OUT  rho,  sound  :  REAL) ; 
END  (DEFINITION)  MODULE  (calcMod). 
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IMPLEMENTATION  MODULE  calcMod; 

FROM  MathMod  IMPORT~POWER,  SIN,  COS,  SORT,  EXP; 

FROM  globalMod  IMPORT  re; 

PROCEDURE  gravCalc  (IN  a:  REAL)  :  REAL; 

BEGIN 

RETURN  32 . 1741*POWER{re/ (a+re)  ,  2.0); 

END  PROCEDURE  {gravCalc}; 

PROCEDURE  densityCalc  (IN  h,  rhoz  :  REAL;  OUT  rho,  sound  :  REAL) ; 
VAR 


t  :  REAL; 


BEGIN 

rho  :=  rhoz  *  EXP (-1 . 0*h/23111 . 0  +  0.294  *  SIN (h/28860 . 0)  +  0.213  * 
SIN(h/86580.0) ) ; 

IF  h  >  0.0 

t  :=  518.688  -  (3 . 56616E-03) *h; 
sound  :=  49.02118  *  SQRT(t); 

IF  h  >  36152.0 

sound  :=  968.08; 

IF  h  >  82345.0 

t  :=  254.988  +  (1 . 64592E-03) *h; 
sound  :=  49.02118  *  SQRT(t); 

IF  h  >  155348.0 

sound  :=  1105.0; 

IF  h  >  262448.0 

sound  :=  846.9; 

IF  h  >  299516.0 

t  :=  -349.812  +  2 . 19456E-03*h 
END  IF; 

END  IF; 

END  IF; 

END  IF; 

END  IF; 

ELSE 

sound  :=  1116.44; 

END  IF; 

END  PROCEDURE  {densityCalc} ; 

END  (IMPLEMENTATION)  MODULE  (calcMod) . 
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Note:  Definition  amd  Implementation  Modules  for  right  and  left  jumpers  are  the  same. 

DEFINITION  MODULE  right JumperMod; 

FROM  lOMod  IMPORT  StreamObj,  FileOSeType (Output) ; 

FROM  globalMod  IMPORT  eType,  matrixType,  vectorType,  encounterType; 

FROM  VortexMod  IMPORT  RightVortexOb j ,  Lef tVortexObj ; 

FROM  VortexControlMod  IMPORT  VortexControl; 

TYPE 


t JumperObj 

=  OBJECT 

fpc, 

g, 

il/ 

gammad, 

12, 

gees , 

lend. 

gmax. 

ipts, 

h. 

ip, 

hmin. 

jl, 

mach. 

j2. 

mass. 

k, 

myTime, 

loop, 

myDrift, 

mpts. 

myDrif tDirection, 

myNumber, 

pb. 

my Plane, 

phi. 

bigloop  : 

INTEGER; 

psi, 

q, 

alt. 

qb. 

alpha. 

qdyn. 

alphad. 

qs. 

beta. 

qsd. 

bxy. 

rad. 

ca. 

rb. 

cao. 

rho. 

caa2. 

rhoz, 

calpha. 

sac. 

cbar. 

sas. 

cbeta. 

sarea. 

cby2. 

salpha. 

cds. 

sbc. 

do. 

sbeta. 

dp. 

sbs. 

cm. 

sinwt. 

cma. 

sound. 

cma2. 

sp. 

cmo. 

phi. 

cmq. 

St. , 

cn. 

t. 

cna. 

theta. 

cna2. 

tar. 

cnb. 

tpoll. 

cnb2. 

tdrift. 

cnq. 

ubl. 

cnr. 

ub2. 

coswt. 

ub3. 

coswtm. 

uel. 

cp. 

ue2. 

cphi. 

ue3. 

csl. 

uxy. 

csn. 

vp. 

ct. 

vpl3. 

cy. 

vpe. 

cyb. 

vpo. 

cyb2. 

w. 

cyr. 

w2. 

dens. 

weight. 

deptime. 

xbod. 

dt. 

xcg. 

dtpr. 

Xdrift, 

dtpoll. 

xlast. 

dtdrift. 

Ydrift  :  REAL 

slength 

:  REAL;  (length  of  suspension  lines) 
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angle  :  REAL;  {the  angle  (in  radians)  which  defines  the  "cone"  of  the 

suspension  lines)  _ 

dcglength  :  REAL;  {distance  from  end  of  suspension  lines  to  paratrooper  c.g.) 

cweight  :  REAL;  {weight  of  canopy) 

sweight  :  REAL;  {weight  of  suspension  lines) 

radius  :  REAL; 
addedmass  :  REAL; 
distcm  :  REAL; 
sysmass  :  REAL; 
paymom  :  REAL; 
distcan  :  REAL; 
distline  :  REAL; 
distpay  :  REAL; 

addDrift  :  BOOLEAN; 

pcds, 

Pt, 

pcdf , 

pm  :  vectorType;  (1X2) 

xe, 

xs, 

xg, 

ue, 

wb, 

vwind, 

temp, 

fb, 

m, 

mb, 

fe, 

uedot, 

wbdot, 

hb  :  vectorType;  (1X3) 

in, 

jn, 

b, 

bn, 

bdot  ;  matrixType;  (3X3) 
lastRightLocation, 

lastLeftLocation  :  eType;  {Dynamic  Array) 

outfile  :  STRING; 

stream  :  StreamObj ; 

encounter  :  encounterType; 

ASK  METHOD  Objinit; 

TELL  METHOD  jump; 

ASK  METHOD  initialize  (IN  sticic  :  INTEGER; 

IN  myPlane  :  INTEGER) ; 

ASK  METHOD  pollVortices  (IN  vortexPlane  :  INTEGER); 

ASK  METHOD  changeDrift; 

ASK  METHOD  findDrift; 

END  OBJECT  { right JumperObj ) ; 

END  (DEFINITION)  MODULE  { right JumperMod) . 
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IMPLEMENTATION  MODULE  right JumperMod; 

FROM  MathMod  ^  IMPORT  EXP,  SIN,  COS,  POWER,  SQRT,  ATAN2 ,  TAN,  pi; 

FROM  VortexMod  IMPORT  RightVortexObj ,  LeftVortexObj ; 

FROM  globalMod  IMPORT  re,  e,  del,  i,  j,  Numberof Planes ,  seedl,  seed2,  seed3, 

repeat; 

FROM  inputMod  IMPORT  streami,  streamE,  streams,  extension,  printTrajectory; 

FROM  globalMod  IMPORT  CrossWindl,  CrossWind2,  CrossWind3,  ShearAltl,  ShearAlt2, 

Headwind,  vsl,  vs2,  vs3; 

FROM  calcMod  IMPORT  gravCalc,  densityCalc; 

FROM  SimMod  IMPORT  SimTime; 

FROM  VortexControlMod  IMPORT  VortexControl,  Airdrop; 

OBJECT  right JumperObj ; 

ASK  METHOD  Objinit; 

BEGIN 

NEW(pcdf,  1..2); 

NEW (pm  ,  1..2); 

NEW ( f b  ,  1 .  .  3 ) ; 

NEW ( f e  ,  1 .  .  3 )  ; 

NEW(hb  ,  1..3); 

NEW (mb  ,  1..3); 

NEW ( temp  ,  1 . . 3 ) ; 

NEW(ue  ,  1..3); 

NEW(uedot,  1..3); 

NEW(vwind,  1..3); 

NEW(wb  ,  1..3); 

NEW(wbdot,  1..3); 

NEW(xe  ,  1..3); 

NEW(xs  ,  1..3); 

NEW(xg  ,  1..3); 


NEW(pcds, 

1. 

.4) 

; 

NEW(pt  , 

1. 

.4) 

/ 

NEW (in  , 

1. 

.3, 

1.  .3); 

NEW ( j  n  , 

1. 

.3, 

1..3); 

NEW(b  , 

1. 

.3, 

1.  .3); 

NEW (bn  , 

1. 

.3, 

1.  .3); 

NEW(bdot, 

1. 

.3, 

1.  .3) ; 

{  -  system  inertial  properties  -  ) 

{  -  parachute-payload  system  weight  (lbs)  =  wight  of  jumper/gear  + 

weight  of  T-IOC  -  ) 

weight  :=  seed3. Normal  (247.0,  24.35); 

mass  :=  weight/32 . 17; 

xcg  :=  2.0;  {  forebody  c.g.  (ft)  in  the  horizontal  } 

xbod  :=  6.0;  {  forebody  length  (ft)  in  the  vertical  } 

FOR  i  :=  1  TO  3 

FOR  j  :=  1  TO  3 

jn[i,j]  :=  0.0; 
in[i,j]  :=  0.0; 

END  (j)  FOR; 

END  (i)  FOR; 

{  -  initial  conditions  -  ) 

FOR  i  :=  1  TO  3 

xe[i]  :=  0.0;  {  (ft)  1:  down  range,  2:  off  range,  3:  altitude  loss 

} 

ue[i]  :=  0.0;  {  (fps)  1:  horizontal  velocity,  2:  lateral  velocity, 

3:  ejection  velocity  positive  down  ) 

wb[i]  :=  0.0;  {  ??? 

} 

vwind[i]  :=  0.0;  {  (fps)  1:  head  (+)  or  tail  (-)  wind,  2:  crosswind, 

3:  ????  ) 

END  FOR; 
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alt  ___ 

hmin 

ue[l] 

ue[31 

theta 

vwind [ 1 ] 

vwincl[2] 

dens 

rhoz 


=  0.0;-  {  altitude  (ft)  } 
=0.0;  {  ground  level  (ft)  } 
=0.0;  {  horizontal  velocity  (fps)  } 
=0.0;  {  ejection  velocity,  positive  down  } 
=0.0;  {  pitch  angle  (deg)  nose  up  positive  ) 
=  0.0;  {  head  (+)  or  tail  (-)  wind  (fps)  ) 
=0.0;  {  crosswind  (fps)  } 
=  0.0;  (  density  (0  for  standard  atms)  in  slug/ft''3  } 
=  0.002378;  {  ????  } 


IF  dens  <>  0.0 

rhoz  :=  dens  *  EXP  (alt/23111.0  -  0.295  *  SIN (alt/28860 . 0)  -  0.213  * 

SIN(alt/86580.0) ) 

END  IF; 


{  -  program  constants  -  } 


) 

} 


dtpr 
dtpoll 
dtdrif t 


every  10  seconds  } 


0.1;  {  print  interval  (sec) 

0.5;  {  poll  vortex  positions  every  0.5  seconds 

5.0;  {  change  drift  angle  +/-  45.0  from  current  drift  angle 


{ 


{ 


{ 


forebody  aerodynamic  coefficients  -  ) 


cbar 

= 

6.0; 

{  reference  length  (ft) 

} 

sarea 

= 

POWER(0 . 5, 2 . 0) *pi;  (  reference  area  (ft"2) 

) 

cna  ^ 

= 

0.0; 

{  normal  force  cn-alpha  (/rad) 

} 

cyb 

= 

0.0; 

{  side  force  cy-beta  (/rad) 

) 

caa2 

= 

0.0; 

{  axial  force  ca-alpha''2  (/rad^2) 

) 

do 

= 

0.0; 

{  roll  torque  coefficient  (dimensionless) 

} 

dp 

0.0; 

{  roll  damping  coefficient  (/rad) 

} 

cma 

= 

-2.0; 

(  pitch  moment  cm-alpha  (/rad) 

) 

cmq 

-200.0; 

{  pitch  damping  (/rad) 

} 

cnb 

= 

0.0; 

{  yaw  moment  cn-beta  (/rad) 

} 

cnr 

0.0; 

{  yaw  damping  (/rad) 

} 

forebody  drag  versus  mach  number  table  -  ) 


mpts  :=  2; 

pm(l]  :=  0.00;  {  mach  number  ) 

pm[21  :=  2.00; 

pcdfdl  :=  0.73+0.06*(360.0  -  weight) /180 . 0;  {  drag  coefficient  } 
pcdf[2I  :=  0.73+0.06*(360.0  -  weight) /180 . 0; 

—  parachute  drag-area  versus  time  table  -  ) 

deptime  :=  0.25;  (  deployment  time  } 
ipts  :=  2; 


IF  deptime  >  0.0 


ipts  :=  ipts  +2;  {  ipts  = 

pt[l] 

=  0.0; 

pt[2] 

=  deptime; 

pt[3] 

=  0.00  +  deptime; 

pt(4] 

=  2.80  +  deptime; 

pcds [1] 

=  0.0; 

pcds [2] 

=  0.0; 

pcds [3] 

=  0.20; 

pcds [ 4 ] 

=  690.0; 

END  IF; 


{  -  convert  EULER  ANGLES  to  direction  cosines  -  ) 


psi 

=  0.0; 

phi 

=  0.0; 

rad 

=  pi/ 180.0 ; 

St. 

:=  SIN  (theta*rad) 

ct 

=  COS  (theta*rad) ; 

sp 

=  SIN  (psi*rad) ; 

cp 

=  COS  (psi*rad); 
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phi; 


phi  :=  SIN  (phi*rad); 
cph^  :=  COS  (phi*racl); 
xe[3]  :=  -1.0  *  alt; 

b[l,ll;=cp*ct; 
b[l,21  :=  sp  *  ct; 
b[l,3]  :=  -1.0  *  St.; 

b[2,l]  :=  -1.0  *  sp  *  cphi  +  cp  *  St.  * 
b[2,2]  :=  cp  *  cphi  +  sp  *  St.  *  phi; 
b[2,3]  :=  ct  *  phi; 

b[3,l]  :=  sp  *  phi  +  cp  *  St.  ♦  cphi; 
b[3,2]  :=  -1.0  *  cp  *  phi  +  sp  *  St.  *  cphi; 
b[3,3]  :=  ct  *  cphi; 

END  {ASK}  METHOD  {Objinit}; 

TELL  METHOD  jump; 

BEGIN 

WHILE  bigloop  =  0 
WAIT  DURATION  dt; 

IF  (-1.0*xe[3])  >  hmin 


SIZE  TO  0.01 

. } 


g  :=  gravCalc{h); 

densityCalc  (h,  rhoz,  rho,  sound); 

IF  cds  =  690.0 

dt  :=  0.01;  {  -  ONCE  CANOPY  INFLATES,  DECREASE  TIME  STEP 

} 

ELSE 

dt  :=  0.001;  {  -  OTHERWISE,  START  WITH  A  SMALLER  STEP  SIZE 

END  IF; 


radius  :=  SQRT(cds/pi) ; 

addedmass  :=  rho* (4 .0/3.0) ♦pi*POWER(radius, 3 . 0) ; 
distcm  := 

(32 . 17*addedmass* (slength*COS (angle) + (4 .0/3 .0) * (radius /pi) tdcglength) 

+ 

cweight* (slength*COS (angle) + (4 . 0/3 . 0) * (radius /pi) +dcglength) 

+  sweight* (0 . 5*slength*COS (angle) tdcglength) )  / 

(32 . 17*addedmass+weight+sweight+cweight) ; 

sysmass  :=  (weight+cweight+sweight) /32 . 17  +  addedmass; 
paymom  :=  ( 1 . 0/12 . 0) *mass* (3 . 0*POWER( (0 . 5*xcg) , 2 . 0)  + 

POWER (xbod, 2.0) ) ; 

distcan  :=  slength*COS (angle)  +  (4 . 0/3 . 0) * (radius/pi)  -  distcm; 
distline  :=  distcm  -  0 . 5*slength*COS (angle) ; 

distpay  :=  distcm; 

in[l,l]  :  = 

(addedmass* { (2. 0/5.0) *POWER (radius, 2 . 0) +POWER(distcan,  2 . 0) )  +  (cweight/32 . 17 ) *POWER{distcan, 
2.0) 

+ 

(sweight/32. 17) *POWER (distline, 2 . 0) +paymom+mass*POWER(distpay, 2.0)) /14 .59; 
in[2,2]  := 

(addedmass* ( (2 . 0/5 . 0) *POWER (radius, 2 . 0) +POWER (distcan, 2 . 0) )+ (cweight/32 . 17) *POWER (distcan, 
2.0) 

+ 

(sweight/32.17) *POWER (distline, 2 .0)+paymom+mass*POWER{distpay, 2.0) ) /14 .59; 

in[3,3]  :=  ( (2 . 0/5 . 0) * (cweight/32 . 17) *POWER (radius, 2 . 0)  + 

(2 .0/3.0) * (rho* (4 . 0/3 . 0) *pi*POWER(radius, 3 . 0) 

*  POWER(radius,2.0) )  +  (0 . 5*mass) ) /14 . 59; 

FOR  i  :=  1  TO  3 

jn[i,i]  :=  1.0  /  in  [i,i]; 

END  FOR; 


IF  (-1.0*xe[3])  <=  ShearAltl 

IF  (-1.0*xe[3])  <=  ShearAlt2 
vwind[2]  :=  vs3; 

ELSE 
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vwind[2]  :=  vs2; 

_  END  IF;- 
ELSE 

vwind[2]  :=  vsl; 

END  IF; 

pb  : =  wb [ 1 ] ; 

qb  :  =  wb  [  2  ]  ; 

rb  :  =  wb  [  3  ]  ; 

uel  :=  ue[l]  -  vwind[l]; 

ue2  ;=  ue[2]  -  vwind[2]; 

ue3  :=  ue[31  -  vwind[3]; 

vp  :=  SQRT(POWER(uel,2.0)  +  P0WER(ue2, 2 . 0)  +  POWER (ue3 , 2 . 0 ) ) 

vpo  :=  vp; 

mach  :=  vp/sound; 

ubl  :=  b[l,l]*uel  +  b(l,21*ue2  +  b[l,3]*ue3; 

ub2  :=  b[2,ll-»uel  +  b[2,2]*ue2  +  b[2,3]*ue3; 

ub3  :=  b[3,l]*uel  +  b(3,21*ue2  +  b[3,3]*ue3; 

vpl3  :=  SQRT(POWER(ubl,2.0)  +  POWER (ub3 , 2 . 0 ) ) ; 

{  -  USE  SIN (ALPHA)  for  ALPHA  and  COS (BETA)  for  BETA  -  ) 

IF  vpo  <  l.OE-06 

vpo  :=  l.OE-06; 

END  IF; 

sbeta  :=  ub2  /  vpo; 
cbeta  :=  vpl3  /  vpo; 
beta  :=  sbeta; 

IF  vpl3  <  l.OE-06 
vpl3  :=  l.OE-06; 

END  IF; 

salpha  :=  ub3  /  vpl3; 
calptia  :=  ubl  /  vpl3; 
alpha  :=  salpha; 

{ - AERODYNAMIC  and  BODY  FORCES  AND  MOMENTS - ) 

{  -  ISOLATED  BODY  AERODYNAMICS  -  ) 

{ - BEGIN  AERO  ROUTINE - ) 

sac  :=  salpha  *  calpha; 

sas  :=  salpha  *  ABS  (salpha); 

sbc  :=  sbeta  *  cbeta; 

sbs  :=  sbeta  *  ABS (sbeta); 

rad  :=  cbar  /  (2.0*vpo); 

cna2  : =  0.0; 

cnq  : =  0.0; 

cyb2  :=  0.0; 

cyr  :=  0.0; 

cmo  : =  0 . 0 ; 

cma2  :=  0.0; 

cnb2  :=  0.0; 

{  -  FOREBODY  AERO-LIFT  DRAG  COEFFIECIENT  -  } 

i  :=  0; 
loop  :=  0; 

WHILE  loop  =  0  (WILL  LOOP  WHEN  mach  >  pm [ip]  UNTIL  i  =  mpts) 
i  :=  i  +  1; 
ip  :=  i  +  1; 

IF  i  =  mpts 

cao  :=  pcdf[2);  (  when  i  =  ipts  ) 
loop  :=  1; 

ELSE 

IF  mach  <=  pm[ip] 

cao  :=  pcdf [i] + (pcdf [ip] -pcdf [i] ) * (mach- 

pm[i] ) / (pm[ip] -pm[i] ) ; 

loop  :=  1; 

END  IF; 
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END  IF; 

END  WHILE; 


cn  :=  cna  *  sac  +  cna2  *  sas  +  cnq  *  qb  *  rad; 

cy  :=  cyb  +  sbc  +  cyb2  *  sbs  +  cyr  *  rb  *  rad; 

ca  :=  cao  +  caa2  *  (1.0  -  POWER(calpha, 2 . 0)  *  POWER (cbeta, 2 . 0 ) ) ; 
csl  :=  do  +  dp  *  pb  *  rad; 

cm  :=  cmo  +  cma  *  sac  +  cma2  *  sas  +  cmq  *  qb  *  rad; 

csn  :=  cnb  *  sbc  +  cnb2  *  sbs  +  cnr  *  rb  *  rad; 

{  -  PARACHUTE  DRAG-AREA  -  } 

cds  : =  0.0; 
i  :  =  0  ; 
loop  :=  0; 

( - THIS  LOOP  INFLATES  THE  PARACHUTE - } 

WHILE  loop  =  0 

IF  t  <  pt[ll 

loop  :=  1; 

ELSE 

i  :=  i  +  1; 

ip  :=  i  +  1; 

IF  i  =  ipts 

cds  :=  pcds(ipts]; 
loop  :=  1; 

ELSE 

IF  t  <=  pt  (ip] 

cds  :=  pcds[il  +  (pcds [ip] -pcds [i] )  *  (t-pt[i])  / 

(pt[ip]-pt[i] ); 

loop  :=  1; 

END  IF; 

END  IF; 

END  IF; 

END  WHILE; 

( - end  aero  ROUTINE - } 

q  :=  0.5  *  rho  *  vp; 
fpc  :=  -1.0  *  q  *  cds; 

qs  :=  0.5  *  rho  *  POWER  (vp,  2.0)  *  sarea; 
qsd  :=  qs  *  cbar; 

fb[l]  :=  -1.0  *  qs  *  ca  +  mass  *  g  *  b[l,3]  +  fpc  *  ubl; 

fb[2]  :=  qs  *  cy  +  mass  *  g  *  b[2,3]  +  fpc  *  ub2; 

fb[3]  :=  -1.0  *  qs  *  cn  +  mass  *  g  *  b[3,3]  +  fpc  *  ub3; 

mb[l]  :=  qsd  *  csl; 

mb[2]  :=  qsd  *  cm  +  fpc  *  ub3  *  (xbod  -  xcg) ; 
mb[3]  :=  qsd  *  csn  -  fpc  *  ub2  *  (xbod  -  xcg); 
gees  :=  -1.0  *  fb[l]  /  (mass  *  g) ; 

IF  ABS  (gees)  >  ABS  (gmax) 
gmax  :=  gees; 

END  IF; 


IF  printTrajectory 

IF  t  >=  tar  {  THEN  PRINT  DATA  } 
tar  :=  tar  +  dtpr; 
h  :=  -1.0  *  xe[3]; 

vpe  :=  SQRT  (  POWER  (uel,2.0)  +  POWER  (ue2,2.0)  +  POWER 

(ue3,2.0)  ); 

qdyn  :=  0.5  *  rho  *  POWER  (vp,2.0); 

bxy  :=  SQRT  (  POWER  (b[l,l],  2.0)  +  POWER  (b[l,2],2.0)  ); 
theta  :=  57.295  *  ATAN2  (  (-1.0  *  b[l,3]),  bxy  ); 
alphad  :=  57.295  *  ATAN2  (  salpha,  calpha  ); 
uxy  :=  SQRT  (  POWER  (ue[l],2.0)  +  POWER  (ue[2],2.0)  ); 
gammad  :=  57.295  *  ATAN2  (  (-1.0  *  ue[3]),  uxy  ); 

{ 


OUTPUT  (myNumber,  "R  ",  SimTime,  "  ”,  h,  "  ”,  xe[l],  "  ", 
xe[2],  "  ”,  vpe,  "  ",  vp,  ”  ”,  mach,  ”  ”,  qdyn,  ”  ”,  gees,  "  ",  gammad,  "  ",  theta,  "  " 
alphad,  "  ",  cds) ; 

) 

ASK  Stream  TO  WriteString  (INTTOSTR(myNumber)  +  "R  "); 
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ASK  stream  TO  WriteString 
^  ASK  stream  TO  WriteString 

ASK  stream  TO  WriteString 
ASK  stream  TO  WriteString 
ASK  stream  TO  WriteString 
ASK  stream  TO  WriteString 
ASK  stream  TO  WriteString 
ASK  stream  TO  WriteString 
ASK  stream  TO  WriteString 
ASK  stream  TO  WriteString 
ASK  stream  TO  WriteString 
ASK  stream  TO  WriteString 
ASK  stream  TO  WriteString 
ASK  stream  TO  WriteLn; 

END  IF; 

END  IF; 


(REALTOSTR(SimTime)  +  "  "); 
(REALTOSTR(h)  +  " 
(REALTOSTR(xe [11  )  +  "  "); 
(REALTOSTR(xe [2] )  +  "  ”); 
(REALTOSTR(vpe)  +  "  "); 
(REALTOSTR(vp)  +  "  "); 
(REALTOSTR(mach)  +  "  "); 
(REALTOSTR(qdyn)  +  "  "); 
(REALTOSTRlgees)  +  "  "); 
(REALTOSTR(gainmad)  +  "  "); 
(REALTOSTR( theta)  +  " 
(REALTOSTR(alphad)  +  "  "); 
(REALTOSTRlcds)  )  ; 


{  EULER  ROTATION  FUNCTION  FOR  DIRECTION  COSINE  PROPOAGATION  } 

w2  :=  POWER(wb(l] ,2.0)  +  POWER (wb [2] , 2 . 0 )  +  POWER (wb [ 3 ] , 2 . 0 ) ; 

w  :=  SQRT  (w2); 

coswt  :=  COS  (w*dt); 

sinwt  :=  SIN  (w*dt) ; 

coswtm  ;=  1.0  -  coswt; 

IF  w2  <  l.OE-12 

w2  :=  l.OE-12; 
w  :=  l.OE-06; 

END  IF; 


{  ANGULAR  MOMENTUM  CROSS  PRODUCT  TERMS  } 

FOR  k  :=  1  TO  3 

hb  [k]  :=  in[k,l]  *  wb[ll  +  in(k,2]  *  wb[21  +  in[k,31  *  wb[3]; 
END  FOR; 

FOR  i  :=  1  TO  3 

11  :=  e(i+l] ; 

12  :=  e [1+2] ; 

temp[i]  :=  wb[il]  *  hb[i2)  -  wb[i2]  *  hb[il]; 

END  FOR; 

(  FORCE  RESOLUTION  TO  EULER  SYSTEM  ) 

{  TRANSLATIONAL  ACCELERATION  AND  DIRECTION  COSINE  ROTATION  } 

FOR  i  :=  1  TO  3 

fe[i]  :=  fb[l)  *  b[l,i]  +  fb(2]  *  b[2,i)  +  fb[3]  *  b[3,i]; 
uedot[i]  :=  fe[i]  /  mass; 

FOR  j  :=  1  TO  3 

bn(i, j]  :=  b[i, j); 
jl  :=  e[j+l]  ; 
j2  :=  e[j+2]; 

bdot[i,j]  ;=  del [i, j ) *coswt  +  wb [i] *wb [ j ] *coswtm/w2  +  ( 
wb [ jl] *d6l [i, j2]  -  wb[ j2] *del [i, jl]  )*  sinwt/w; 

END  FOR; 

{  ANGULAR  ACCELERATION  IN  BODY  AXES  ) 

wbdot[i]  :=  jn(i, 1] * (mb[l] -temp[l] )  +  jn [i, 2] * (mb [2] -temp  [2] )  + 

jn[i,3]*(mb[3]  -temp  [3]  )  ; 

END  FOR; 


( - INTEGRALS - ) 

t  :=  t  +  dt; 

{  t  :=  SimTime  -  myTime;} 

xlast  :=  xe [ 1] ; 

FOR  i  :=  1  TO  3 
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xe[il  :=  xe[i]  +  dt* (ue [i] +0 . 5*dt*uedot [i] ) ; 

„  ue[il  ue(i]  +  dt*uedot(i]; 
wb[i]  :=  wb[i]  +  dt*wbdot[i]; 

FOR  j  ;=  1  TO  3 

b[i,j]  :=  bdot[i,l]*bn[l, j]  +  bdot [i, 2] *bn[2, j ]  + 

bdot [i, 3 ] *bn [ 3 , j ] ; 

END  FOR; 

END  FOR; 

(  -  INDUCE  A  DRIFT  DIRECTION  AND  VELOCITY  ON  THE  PARATROOP 


IF  t  >=  tdrift 

tdrift  :=  t  +  dtdrift; 

ASK  SELF  TO  changeDrift; 

END  IF; 

IF  t  >=  6.5 

xe[l]  :=  xe[l)  +  Xdrift*dt; 
xet2]  ;=  xe[2]  +  Ydrift*dt; 

END  IF; 

IF  t  <  4.1 

xe[2]  :=  xe(2]  -  9.25/4100.0; 

END  IF; 

{  -  UPDATING  MOVING  AND  GROUND  COORDINATE  SYSTEMS  -  } 

xs[l]  :=  Airdrop. Information [1] .xg  -  xe[l]; 
xs[2]  :=  xe[21  ; 
xs [ 3 ]  ; =  xe [ 3 ] ; 

xg[l]  :=  Airdrop. InformationCl] .xg  -  xs[l]; 
xg[2]  :=  xe(2]; 
xg [ 3 1  : =  xe [ 3 ] ; 

{ - POLL  ALL  VORTICES  FOR  MISSED  DISTANCE - } 

IF  t  >=  tpoll 

IF  cds  >=  pods [4] 

tpoll  :=  tpoll  +  dtpoll; 

FOR  i  :=  1  TO  myPlane-1; 

ASK  SELF  TO  pollVortices  (i); 

END  FOR; 

END  IF; 

END  IF; 

ELSE  { - when  -xe[3]  <=  hmin,  THEN  PRINT  DATA  FOR  LAST  TIME - ) 

h  :=  -1.0  *  xe[3]; 

vpe  ;=  SORT  (  POWER  {uel,2.0)  +  POWER  (ue2,2.0)  +  POWER  (ue3,2.0) 

)  ; 

qdyn  :=  0.5  *  rho  *  POWER  (vp,2.0); 

bxy  :=  SQRT  (  POWER  (b(l,l),2.0)  +  POWER  (b[l,2],2.0)  ); 

theta  :=  57.295  *  ATAN2  (  (-1.0  *  b[l,3]),  bxy  ); 

alphad  :=  57.295  *  ATAN2  (  salpha,  calpha  ); 

uxy  :=  SQRT  (  POWER(ue[l] ,2.0)  +  POWER (ue [2] , 2 . 0)  ); 

gammad  :=  57.295  *  ATAN2  (  (-1.0  *  ue[3]),  uxy  ); 

OUTPUT  (myNumber,  "R  ",  myPlane,  "  ",  SimTime,  ”  ",  h,  "  ",  xe[l],  " 
",  xe[2],  "  ",  vpe,  "  ",  vp,  "  ",  mach,  "  ",  qdyn,  "  ",  gees,  "  ",  gammad,  "  ",  theta,  " 
",  alphad,  "  ",  cds); 

ASK  Streams  TO  WriteString  (INTTOSTR( repeat)  +  "  "); 

ASK  Streams  TO  WriteString  (INTTOSTR (myNumber)  +  "R  ") ; 

ASK  Streams  TO  WriteString  ( INTTOSTR (myPlane)  +  "  "); 

ASK  Streams  TO  WriteString  (REALTOSTR( SimTime)  +  "  "); 

ASK  Streams  TO  WriteString  (REALTOSTR (h)  +  "  "); 

ASK  Streams  TO  WriteString  (REALTOSTR (xe [1] )  +  "  "); 

ASK  Streams  TO  WriteString  (REALTOSTR (xe [2] )  +  "  "); 

ASK  Streams  TO  WriteString  (REALTOSTR (vpe)  +  "  ") ; 

ASK  Streams  TO  WriteString  (REALTOSTR (vp)  +  "  ") ; 

ASK  Streams  TO  WriteString  (REALTOSTR (mach)  +  "  "); 

ASK  Streams  TO  WriteString  (REALTOSTR (qdyn)  +  "  "); 
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ASK 

Streams 

TO 

ASK 

Streams 

TO 

ASK 

streams 

TO 

ASK 

streams 

TO 

ASK 

streams 

TO 

ASK 

Streams 

TO 

WriteString 

Writestring 

WriteString 

WriteString 

WriteString 

WriteLn; 


(REALTOSTR(gees) 

( REALTOSTR  ( gainmad ) 
(REALTOSTRI theta) 
(REALTOSTR(alphad) 
( REALTOSTR (cds)  )  ; 


+ 

+ 

+ 


+ 


IF  printTrajectory 


)  ; 
)  ; 
)  ; 
)  ; 


ASK  stream  TO  WriteString  ( INTTOSTR (myNumber )  +  "R  "); 
ASK  stream  TO  WriteString  { REALTOSTR (SimTime)  +  "  "); 


ASK  stream  TO  WriteString  (REALTOSTR (h)  +  "  ") 
ASK  stream  TO  WriteString  {REALTOSTR(xe [1] )  +  "  ") 
ASK  stream  TO  WriteString  (REALTOSTR (xe [2] )  +  "  ") 
ASK  stream  TO  WriteString  ( REALTOSTR (vpe)  +  "  ") 
ASK  stream  TO  WriteString  ( REALTOSTR (vp)  +  "  ") 
ASK  stream  TO  WriteString  (REALTOSTR (mach)  +  "  ") 
ASK  stream  TO  WriteString  ( REALTOSTR (qdyn)  +  "  ") 
ASK  stream  TO  WriteString  (REALTOSTR(gees)  +  "  ") 
ASK  Stream  TO  WriteString  ( REALTOSTR (gammad)  +  "  ") 
ASK  Stream  TO  WriteString  (REALTOSTR(theta)  +  "  ") 
ASK  stream  TO  WriteString  (REALTOSTR (alphad)  +  "  ") 
ASK  stream  TO  WriteString  (REALTOSTR (cds) ) ; 

ASK  stream  TO  WriteLn; 


END  IF; 

bigloop  :=  1; 

END  IF  (hmin); 

END  WAIT; 

END  WHILE; 

IF  printTrajectory 

ASK  stream  TO  Close; 
DISPOSE  (stream); 

END  IF; 


DISPOSE  (SELF); 


END  {ASK}  METHOD  (jump); 


ASK  METHOD  initialize  (IN  stick  :  INTEGER; 

IN  Counter  :  INTEGER)  ; 


BEGIN 


myPlane 

myNumber 

gees 

cds 

myTime 

t 

gmax 

tar 

tpoll 

tdrift 

alt 

h 

bigloop 

myDrift 

myDriftDir action 

xe  [1] 

xlast 

xe[2] 

alt 

xe  [3] 

ue  [1] 


:=  Counter; 

:=  stick; 

:=  0.0; 

:=  0.0; 

:=  SimTime; 

:=  SimTime-myTime; 

:=  0.0; 

:=  t; 

:=  t; 

:=  t  +  dtdrift; 

:=  -1.0  *  xe[3]; 

:=  alt; 

:=  0; 

:=  seedl .UniformReal  (0.0,  4.0); 

:=  seed2 .UniformReal  (0.0,  360.0); 

:=  Airdrop. Information [myPlane] .xg; 

: =  xe [ 1 ] ; 

:=  Airdrop. Information [myPlane] . yg  +  9.25; 
:=  Airdrop. Information [myPlane] .altitude; 
:=  -1.0  *  alt; 

Airdrop. Information [myPlane] .vf; 


vwind[l] 


:=  Headwind;  [FROM  globalMod) 
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NEW  (lastRightLocation,  1. .myPlane-l) ; 
NEW  ^lastLeftLocation,  1 . .myPlane-1)  ; 


FOR  i  :=  1  TO  myPlane-1 

lastRightLocation[i)  :=  1; 
lastLeftLocation [i]  :=  1; 

END  FOR; 

ASK  SELF  TO  findDrift; 

IF  printTrajectory 
NEW  (stream); 

outfile  :=  "RJ"  +  INTTOSTR(myPlane)  +  INTTOSTR{myNumber)  +  extension  + 

.mat"  ; 

ASK  Stream  TO  Open  (outfile.  Output) ; 

END  IF; 

ASK  streami  TO  WriteString  (INTTOSTR( repeat)  +  "  "); 

ASK  streami  TO  WriteString  ( INTTOSTR (myPlane )  +  "  ") ; 

ASK  streami  TO  WriteString  (REALTOSTR(SimTime)  +  "  "); 

ASK  streami  TO  WriteString  (INTTOSTR(stick)  +  "R  "); 

ASK  streami  TO  WriteString  (REALTOSTR(weight)  +  ”  "); 

ASK  streami  TO  WriteString  (REALTOSTR (xe [ 1 ] )  +  "  "); 

ASK  streami  TO  WriteString  (REALTOSTR (xe [2] )  +  "  "); 

ASK  streami  TO  WriteString  (REALTOSTR (alt)  +  ”  "); 

ASK  streami  TO  WriteLn; 

END  METHOD  (initialize); 

ASK  METHOD  pollVortices  (IN  vortexPlane  ;  INTEGER); 

VAR 


X,  xcordl,  xcord2, 

vvx,  vjx 

REAL; 

y,  ycordl,  ycord2, 

vvy,  vjy 

REAL; 

z,  zcordl,  zcord2, 

wz,  viz 

REAL; 

vjdistance,  vvdistance,  distance 

REAL; 

projection 

REAL; 

i,  location 

INTEGER; 

check 

BOOLEAN; 

StartRightSearch 

INTEGER; 

startLeftSearch 

INTEGER; 

BEGIN 

checlc  :=  FALSE; 
location  ;=  0; 

startRightSearch  ;=  lastRightLocation [vortexPlane] ; 

{ - POLL  RIGHT  VORTEX - ) 

FOR  i  :=  StartRightSearch  TO  Airdrop. Information [vortexPlane] .NumberOfSteps 
IF  location  =  0 
IF 

(ABS (Airdrop. Information [vortexPlane] .C17.RightVortex.CompletePosition[i] .xCord-xs[l] ) )  <= 
50.0 

lastRightLocation [vortexPlane]  :=  i; 

IF 

(ABS (Airdrop . Information [vortexPlane] .017 .RightVortex. Complete Posi tion [i] . zCord+xs [3] ) )  <= 
50.0 

lastRightLocation [vortexPlane]  :=  i; 

IF 

(ABS (Airdrop . Information [vortexPlane] .017 . RightVortex.CompletePosition [i] . yCord-xs [2] ) )  <= 
50.0 

check  :=  TRUE; 

location  :=  i; 

IF  location  = 

Airdrop . Information [vortexPlane] .NumberOfSteps 

check  :=  FALSE; 

END  IF; 

lastRightLocation[vortexPlane]  :=  i; 

ELSE 
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location  :=  i; 

_  -  END  IF; 

~  ELSE 

location  :=  i; 

END  IF; 

END  IF; 

ELSE 

EXIT; 

END  IF; 

END  FOR; 

IF  check 

xcordl  := 

Airdrop .  Information  [vortex Plane]  .  C17  .  RightVortex.CoinpletePosition  [location]  .xCord; 
ycordl  := 

Airdrop . Information [vortexPlane] .C17 . RightVortex .Complete Posit ion [location] . yCord; 
zcordl  := 

Airdrop . Information [vortexPlane] .C17 . RightVortex. Complete Posit ion [location] . zCord; 
xcord2  := 

Airdrop . Information [vortexPlane] . Cl 7 .RightVortex. CompletePos it ion [location+1] . xCord; 
ycord2  := 

Airdrop . Information [vortexPlane] .C17 . RightVortex .CompletePosition [location+1 ] . yCord; 
zcord2  := 

Airdrop . Information [vortexPlane] . C17 .RightVortex. Complete Posit ion [location+1] . zCord; 

wx  :=  xcord2 -xcordl; 
wy  :=  ycord2 -ycordl; 
wz  :=  zcord2-zcordl; 

vjx  :=  xs [1] -xcordl; 
vjy  ;=  xs [2] -ycordl; 
vjz  :=  -1 . 0*xs [3] -zcordl; 

vjdistance  :=  SQRT (POWER(vjx, 2 . 0) +POWER (vjy, 2 . 0) +POWER (vj z, 2 . 0) ) ; 
wdistance  :=  SQRT(POWER(wx,2.0)+POWER{vvy,2.0)+POWER(wz,2.0)  )  ; 
projection  :=  (vjx*vvx  +  vjy*wy  +  vjz*vvz) /wdistance; 
distance  :=  SQRT(POWER(vjdistance, 2. 0) -POWER (projection, 2 . 0) ) ; 

IF  distance  <= 

MAXOF (Airdrop. Information [vortexPlane] . C17. RightVortex. CompletePosi tion[ location] .radius, 
Airdrop. Information [vortexPlane] .C17 .RightVortex .CompletePosition [location+1] .radius) 

OUTPUT  (myNumber,  "R  ”,  myPlane,  ”  RV  ",  vortexPlane,  "  ”,  - 
1.0*xe[31,  ”  ",  distance,  ”  ",  location,  "  ",  SimTime) ; 

ASK  streamE  TO  WriteString  (INTTOSTR(repeat)  +  "  ") ; 

ASK  StreamE  TO  WriteString  (INTTOSTR(myNumber)  +  "R  ") ; 

ASK  StreamE  TO  WriteString  (INTTOSTR(myPlane)  +  "  RV  "); 

ASK  StreamE  TO  WriteString  (INTTOSTR (vortexPlane)  +  "  "); 

ASK  StreamE  TO  WriteString  (REALTOSTR (-1 . 0*xe [3] )  +  "  "); 

ASK  StreamE  TO  WriteString  (REALTOSTR (distance)  +  "  "); 

ASK  StreamE  TO  WriteString  (REALTOSTR(location)  +  "  ") ; 

ASK  StreamE  TO  WriteString  (REALTOSTR(SimTime) ) ; 

ASK  StreamE  TO  WriteLn; 

END  IF; 

END  IF; 

check  :=  FALSE; 
location  :=  0; 

startLeftSearch  :=  lastLeftLocation [vortexPlane] ; 

{ - POLL  LEFT  VORTEX - } 

FOR  i  :=  StartLeftSearch  TO  Airdrop. Information [vortexPlane] .NumberOf Steps 
IF  location  =  0  ' 

IF 

(ABS (Airdrop . Information [vortexPlane] .C17 .LeftVortex. CompletePosition [i] .xCord-xs [1] ) )  <= 
50.0 

lastLeftLocation[vortexPlane]  :=  i; 


IF 

(ABS (Airdrop . Information [vortexPiane] .C17 .LeftVortex .Complete Posit ion [i] . zCord+xs [ 3] ) )  < 
50.0 

lastRightLocation [vortexPiane]  :=  i; 

IF 

(ABS (Airdrop . Information [vortexPiane] .C17 . LeftVortex.CompletePosition [i] .yCord-xs [2] ) )  < 
50.0 

check  :=  TRUE; 

location  :=  i; 

IF  location  = 

Airdrop . Information [vortexPiane] .NumberOf Steps 

check  :=  FALSE; 

END  IF; 

lastRightLocation[vortexPlane]  :=  i; 

ELSE 

location  :=  i; 

END  IF; 

ELSE 

location  :=  i; 

END  IF; 

END  IF; 

ELSE 

EXIT; 

END  IF; 

END  FOR; 

IF  check 

xcordl  := 

Airdrop . Information [vortexPiane] .C17 . LeftVortex.CompletePosition [location] .xCord; 
ycordl  := 

Airdrop . Information [vortexPiane] . C17 . LeftVortex.CompletePosition [location] . yCord; 
zcordl  := 

Airdrop . Information [vortexPiane] .C17 .LeftVortex.CompletePosition [location] . zCord; 
xcord2  ;= 

Airdrop . Information [vortexPiane] . C17 . LeftVortex .CompletePos it ion [location+1]  . xCord; 
ycord2  := 

Airdrop. Information [vortexPiane] .C17 .LeftVortex.CompletePosition [location+1] .yCord; 
zcord2  ;= 

Airdrop . Information [vortexPiane] . C17 . LeftVortex . Complete Position [ location+1 ] . zCord; 

wx  :=  xcord2-xcordl; 
wy  :=  ycord2-ycordl; 
wz  :=  zcord2-zcordl; 

vjx  :=  xs [1] -xcordl; 
vjy  :=  xs [2] -ycordl; 
vj z  :=  -1 . 0*xs [3] -zcordl; 

vjdistance  :=  SQRT (POWER(vjx, 2 . 0) +POWER(vjy, 2 . 0) +POWER(vj z, 2 . 0) ) ; 
wdistance  :=  SQRT  (POWER(vvx,  2 . 0) +POWER  ( wy,  2 . 0) +POWER(vvz,  2 . 0)  )  ; 
projection  :=  (vjx*wx  +  vjy*wy  +  vjz*wz) /wdistance; 
distance  :=  SQRT (POWER(vjdistance, 2 . 0) -POWER(projection, 2 . 0) ) ; 

IF  distance  <= 

MAXOF( Airdrop. Information [vortexPiane] .C17. LeftVortex.CompletePosition [location] .radius, 
Airdrop . Information [vortexPiane] .C17 . LeftVortex.CompletePosition [location+1] . radius ) 

OUTPUT  (rayNumber,  "R  ",  myPlane,  "  LV  ",  vortexPiane,  "  ",  - 
1.0*xe[3],  "  ",  distance,  "  ",  location,  "  ",  SimTime) ; 

ASK  streamE  TO  WriteString  (INTTOSTR( repeat)  +  "  ")  ; 

ASK  StreamE  TO  WriteString  {INTTOSTR (myNumber )  +  "R  ") ; 

ASK  StreamE  TO  WriteString  ( INTTOSTR (myPlane)  +  "  LV  ") ; 

ASK  StreamE  TO  WriteString  (INTTOSTR(vortexPlane)  +  "  ") ; 

ASK  StreamE  TO  WriteString  (REALTOSTR(-l . 0*xe [3] )  +  "  "); 

ASK  StreamE  TO  WriteString  (REALTOSTR(distance)  +  "  "); 

ASK  StreamE  TO  WriteString  (REALTOSTR(location)  +  "  "); 

ASK  StreamE  TO  WriteString  ( REALTOSTR( SimTime )) ; 

ASK  StreamE  TO  WriteLn; 

END  IF; 

END  IF; 
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END  METHOD  {pollVortices ) ; 
ASK  METHOD  changeDrift; 

BEGIN 


myDrift  :=  seedl.UniformReal  (0.0,  4.0); 

myDriftDirection  :=  seed2. Normal  (myDriftDirection,  2.8125) 
ASK  SELF  TO  findDrift; 

END  {ASK}  METHOD  (changeDrift); 

ASK  METHOD  findDrift; 

BEGIN 


Xdrift  :=  myDrift  *  COS  (myDriftDirection*pi/180 . 0 ) ; 
Ydrift  ;=  myDrift  *  SIN  (myDriftDirection*pi/180 . 0) ; 

END  (ASK)  METHOD  {findDrift}; 

END  OBJECT  { right JumperObj } ; 

END  (IMPLEMENTATION)  MODULE  ( right JumperMod) . 
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APPENDIX  E 

Propagation  Time  Step  Comparisons 
(No  Winds) 


Figure  31.  Altitude  Vs.  Time 
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Figure  33.  Altitude  Vs.  Down  Range 
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Figure  36.  Mach  Number  Vs.  Time 


AXIAL  ACCELERATION  DYNAMIC  PRESSURE 


Figure  37.  Dynamic  Pressure  Vs.  Time 


Figure  38.  Axial  Acceleration  Vs.  Time 
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Figure  42.  Altitude  Vs.  Time 
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Figure  46.  Descent  Velocity  Vs.  Time 
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Figure  48.  Mach  Number  Vs.  Time 


Figure  49.  Dynamic  Pressure  Vs.  Time 
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Figure  50.  Axial  Acceleration  Vs.  Time 


Figure  51.  Trajectory  Angle  Vs.  Time 
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TIME 


Figure  53.  System  Rotation  Vs.  Time 
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APPENDIX  F 

Experimental  Design  Analysis  Results 


Table  14 

Design  Point  Description 


Design  Point 

Seed  Used  Trail  Distance  (ft) 

Real  Coded 

Lateral  Separation  (ft) 

Real  Coded 

1 

1 

- 

0 

- 

2 

2 

* 

500 

+ 

3 

3 

+ 

0 

- 

4 

4 

32,000 

+ 

500 

+ 

5 

5 

23,500 

0 

250 

0 

6 

6 

15,000 

- 

250 

0 

7 

7 

23,500 

0 

0 

- 

8 

8 

32,000 

+ 

250 

0 

9 

9 

23,500 

0 

500 

+ 
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Table  15 
Design  Point  1 

Mean  Distance  From  Assembly  Areas 
(Standard  Deviation  Italicized) 


Run 

1 

Company 

2 

3 

1 

430.0703 

66.9743 

360.9418 

35.0532 

389.4421 

56.0624 

2 

435.738! 

76. 7285 

359.4109 

33.2316 

384.7579 

54.573 

3 

437.2077 

73.434 

361.274 

36.7194 

384.6344 

54.731 

4 

439.495 

76.6081 

362.7578 

34.8662 

389.312 

55.056 

5 

432.4826 

78.5673 

362.8335 

37.6427 

387.8803 

57.4444 

6 

440.9464 

74.4567 

364.7747 

33.7632 

389.0972 

55.2675 

7 

435.1913 

71.6705 

360.851 

37.2194 

385.1381 

59.8104 

8 

436.1233 

73.086 

361.5649 

32.6076 

389.3502 

59.514 

9 

439.5211 

74.6969 

364.8641 

34.6763 

387.8497 

62.3549 

10 

436.3754 

70.6358 

366.1484 

33.2577 

393.4379 

60.258 

11 

433.4462 

75.7921 

366.547 

26.8489 

388.3404 

56.75-^7 

12 

438.1389 

76.0421 

362.4927 

38.9077 

384.8852 

61.7537 
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Table  16 
Design  Point  2 

Mean  Distance  From  Assembly  Areas 
(Standard  Deviation  Italicized) 


Run 

1 

Company 

2 

3 

1 

503.6917 

95.4484 

291.8136 

83.1988 

320.3443 

99.9413 

2 

505.4644 

1 00. 306 1 

286.7763 

84.0816 

320.9853 

99.535 

3 

504.9273 

97.159 

289.429 

85.1705 

313.3477 

96.2883 

4 

506.7988 

102.3825 

285.4877 

87.3584 

320.8816 

94.9966 

5 

502.9726 

94.5522 

280.5451 

80.9113 

324.3242 

100.4835 

6 

504.9503 

95.2559 

292.0325 

82.4191 

319.6078 

98.1254 

7 

501.9642 

100.7286 

286.1554 

86.4196 

314.252 

98.6861 

8 

502.0507 

102.2936 

287.9634 

82.949 

320.9301 

92.7337 

9 

504.6824 

97.807 

288.4857 

85.9369 

320.0993 

95.3397 

10 

507.5064 

96.7031 

290.7363 

5J.60P7 

316.6711 

94.6602 

11 

507.2636 

97.0514 

283.3885 

81.6215 

321.9584 

97.572/ 

12 

507.2985 

99.1181 

284.5784 

82.6224 

317.1065 

90.1595 

F-3 


Table  17 
Design  Point  3 

Mean  Distance  From  Assembly  Areas 
(Standard  Deviation  Italicized) 


Run 

1 

Company 

2 

3 

1 

430.4193 

76.6798 

362.1975 

33.6641 

395.0287 

63.3109 

2 

429.7015 

70.1104 

364.5905 

32.5386 

391.5208 

57.8174 

3 

428.6318 

75.4378 

358.3023 

33.2877 

395.2355 

64.0831 

4 

428.983 

74.1196 

361.1818 

29.0871 

397.8057 

61.7333 

5 

428.9359 

70.7156 

358.9821 

28.5781 

391.6958 

60.588 

6 

426.9134 

74.828 

364.7675 

31.3469 

389.6335 

64.0662 

7 

427.2351 

67.914 

361.5911 

29.0319 

390.125 

59.5754 

8 

424.3669 

75.1556 

359.0497 

i2.P695 

390.1488 

58.6665 

9 

430.1527 

69.624 

365.4308 

26.3724 

389.9997 

63.5268 

10 

423.7409 

73.8034 

360.0008 

28.0813 

394.0922 

(55.P(553 

11 

428.9483 

74.7155 

358.8766 

30.6835 

394.8331 

5<5.7i7P 

12 

428.4258 

72.8381 

361.8863 

32.1484 

393.07 

63.8151 
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Table  18 
Design  Point  4 

Mean  Distance  From  Assembly  Areas 
(Standard  Deviation  Italicized) 


Run 

1 

Company 

2 

3 

1 

499.6217 
94. 7887 

284.3822 

81.0878 

328.0279 

106.4739 

2 

496.7935 
99. 7158 

287.3618 

82.0323 

330.9178 

103.4148 

3 

494.2892 

98.6101 

287.2261 

86.2168 

328.161 

102.6164 

4 

495.6269 

95.4149 

283.0487 

84.8408 

329.2329 

98.2828 

5 

497.7153 

98. 7909 

287.8541 

83.7419 

324.6052 

100.9202 

6 

498.8645 

94.7652 

285.3582 

M7337 

326.361 

97.9557 

7 

497.6243 

93.7862 

287.7428 

86.3785 

324.3128 

102.5909 

8 

502.3004 

95.6655 

285.5206 

87.4673 

326.4205 

95.9542 

9 

498.8575 

95.2053 

284.9746 

<55.7642 

326.3337 

101.6063 

10 

496.2071 

93.6687 

285.2436 

84.3306 

327.2455 

106.2513 

11 

491.7247 

94.6025 

284.9658 

81.9371 

329.3919 

104.9135 

12 

496.7456 

99.7924 

283.941 

83.141 

328.4583 

97.1479 
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Table  19 
Design  Point  5 

Mean  Distance  From  Assembly  Areas 
(Standard  Deviation  Italicized) 


Run 

1 

Company 

2 

3 

1 

468.1302 

79.034 

319.4243 

53.2432 

354.674 

72.7081 

2 

502.6151 

63.0635 

283.7802 

43.202 

318.1862 

77.1984 

3 

467.n&2 

80.9138 

321.6672 

50.5145 

354.2458 

70.0561 

4 

461.9555 

80.5971 

322.3655 

50.6149 

360.5623 

73.8028 

5 

463.5697 

77.5801 

321.7715 

48.3536 

356.1365 

73.9633 

6 

465.8508 

74.9795 

323.945 

48.0725 

360.0028 

75.4596 

7 

465.0878 

77.3808 

322.6825 

51.0495 

8 

470.5861 

79.7612 

324.7836 

53.0919 

352.4648 

71.1524 

9 

467.7505 

78.0351 

320.5336 

54.3631 

347.3068 

74.4806 

10 

464.2716 

79.9024 

■325.ini 

57.1777 

355.7873 

69.5746 

11 

462.4232 

80.4521 

325.5173 

51.8559 

357.6986 

71.1894 

12 

462.7124 

79.1749 

323.9769 

51.9871 

352.3883 

70.2033 
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Table  20 
Design  Point  6 

Mean  Distance  From  Assembly  Areas 
(Standard  Deviation  Italicized) 


Run 

1 

Company 

2 

3 

1 

466.2106 

82.099 

321.3344 

55.28 

354.1209 

73,5537 

2 

465.632 

78.3881 

322.7295 

53.9951 

349.6821 

76.4318 

3 

471.0641 

76.9606 

325.9871 

53.4872 

351.8965 

73.  P2 

4 

464.0373 

80.1052 

322.2092 

52.102 

356.4903 

69.5845 

5 

469.5666 

80.8256 

324.6656 

51.5468 

351.3001 

76.5468 

6 

468.87 

78.3476 

326.2134 

52.0792 

352.3676 

75.2512 

7 

469.3715 

78.6952 

319.1993 

51.3488 

351.8846 

74.4209 

8 

462.1089 

74.1749 

323.1526 

50.8189 

352.1188 

77.4046 

9 

468.4109 

81.6074 

322.657 

52.308 

347.4057 

70.4199 

10 

475.9776 

81.2664 

324.2565 

52.4166 

357.0102 

75.7948 

11 

468.8593 

81.1118 

324.9368 

51.1617 

353.1134 

70.5241 

12 

468.9888 

77.6527 

324.6072 

52.3715 

349.1569 

65.3615 
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Table  21 
Design  Point  7 

Mean  Distance  From  Assembly  Areas 
(Standard  Deviation  Italicized) 


Run 

1 

Company 

2 

3 

1 

434.2541 

82.2195 

358.2854 

25.0726 

390.6165 

59.428 

2 

431.828 

74. 7502 

361.756 

31.3976 

392.561 

52.9896 

3 

430.5695 

75.9724 

355.67 

32.4983 

388.312 

58.6882 

4 

429.3265 

75.4169 

358.9588 

28.8031 

392.4118 

60.5126 

5 

430.1651 

68.8505 

363.8356 

32.6825 

391.7858 

59.0724 

6 

436.1033 

68.4442 

357.7043 

34.3462 

392.7163 

59.0649 

7 

434.2491 

80.3527 

362.9033 

34.535 

390.6143 

58.0623 

8 

431.41 

72.5305 

357.7596 

33.9889 

391.1626 

57.4189 

9 

433.5167 

73.5515 

356.8865 

30.6631 

397.0389 

57.9238 

10 

434.4781 
73.1 10 1 

362.4084 

31.4961 

389.9141 

60.8326 

11 

434.967 

77.0348 

361.259 

32.9094 

394.661 

57.6345 

12 

435.7988 

73.5451 

360.1292 

27.3131 

388.6504 

57.2556 
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Table  22 
Design  Point  8 

Mean  Distance  From  Assembly  Areas 
(Standard  Deviation  Italicized) 


Run 

1 

Company 

2 

3 

1 

459.4309 

72.5859 

319.2125 

50.0239 

355.7095 

7(5.9P29 

2 

459.2828 
78. 1941 

325.5855 

52.9487 

360.2764 

77.9914 

3 

465.1589 

80.4955 

319.3877 

51.5492 

362.3887 

76.8498 

4 

458.7631 

75.5444 

321.2724 

51.9716 

359.7524 

76.7526 

5 

459.2312 

76.3075 

322.9782 

48.9529 

357.3255 

74.8719 

6 

461.3052 

74.6898 

321.6286 

52.5543 

365.8863 

76.4414 

7 

457.8445 

81.1468 

324.0579 

51.3617 

358.9947 

79.3221 

8 

457.9037 

79.7287 

322.0138 

48.0619 

3S9.A9%1 

75.4262 

9 

466.3348 

78.7731 

326.1008 

49.9689 

360.3373 

69.7964 

10 

461.259 

79.2548 

320.7537 

46.7226 

361.4708 

74.9443 

11 

463.7841 

76.6226 

321.9678 

50.8494 

359.117 

77.2573 

12 

461.5772 

80.1016 

322.0916 

49.8423 

361.9966 

79.7121 
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Table  23 
Design  Point  9 

Mean  Distance  From  Assembly  Areas 
(Standard  Deviation  Italicized) 


Run 

1 

Company 

2 

3 

1 

501.7365 
99. 7392 

290.1972 

83.4302 

324.4814 

100.7729 

2 

499.0175 

94.1012 

284.8344 

80.6012 

321.6855 

102.6073 

3 

500.6175. 

99.8688 

290.4487 

82.3488 

321.9107 

101.3741 

4 

499.6213 

97.2009 

286.9653 

85.7468 

324.1818 

P7. 7752 

5 

492.5589 

1 01. 320 5 

291.6181 

88.3113 

318.7303 

100.0347 

6 

498.8291 

95.0091 

282.1408 

88.7394 

319.6981 

93.9041 

7 

504.2666 

96.6847 

282.5252 

83.585 

320.3021 

101.5246 

8 

505.8378 

97.7124 

290.7543 

84.4505 

319.686 

99.525 

9 

500.845 

98.2773 

288.8181 

85.3541 

328.0465 

101.6197 

10 

496.5483 

1 00. 1 809 

283.4779 

86.3349 

322.9774 

98.4959 

11 

498.5302 

99.3958 

286.4457 

87.6416 

324.5565 

98.5074 

12 

501.3392 

95.3014 

284.8207 

86.069 

322.2554 

106.6055 
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Table  24 

Summary  of  Results  For 
Mean  Distance  to  Assembly  Areas 


Design  Point 

1 

Company 

2 

3 

1 

Overall  Mean 

436.228 

3.154 

362.8717 

2.2571 

387.8438 

2.6283 

Mean  STD 

74.0577 

3.1681 

34.5662 

3.1314 

57.795 

2.8214 

2 

Overall  Mean 

504.9642 

1.9993 

287.2827 

3.4923 

319.209 

3.2287 

(-,+) 

Mean  STD 

98.2338 
2.70 1 1 

83.8582 

1.9889 

96.8437 

3.1652 

3 

Overall  Mean 

428.0379 

2.1276 

361.4047 

2.4866 

392.7657 

2.641 

(+,-) 

Mean  STD 

72.9951 

2.7497 

30.6491 

2.3736 

61.6572 

2.442 

4 

Overall  Mean 

497.1976 

2.6989 

285.635 

1.5706 

327.4557 

1.961 

(+5+) 

Mean  STD 

96.2339 

2.3061 

84.306 

2.0721 

101.5088 

3.5443 

5 

Overall  Mean 

468.5109 

11.0566 

319.6805 

11.4708 

351.9928 

11.2221 

(0,0) 

Mean  STD 

77.5729 

4.8712 

51.1272 

3.5312 

72.7923 

2.3825 

6 

Overall  Mean 

468.2581 

3.5474 

323.4957 

2.0236 

352.2123 

2.791 

(-,0) 

Mean  STD 

79.2695 

2.3322 

52.4096 

1.2823 

73.2678 

3.5897 

7 

Overall  Mean 

433.0555 

2.3006 

359.7963 

2.6265 

391.7037 

2.4526 

(0,-) 

Mean  STD 

74.6482 

4.0372 

31.3088 

2.P275 

58.2386 

2.014 

8 

Overall  Mean 

460.9896 

2.8137 

322.2542 

2.1514 

360.2295 

2.555P 

(+,0) 

Mean  STD 

77.7871 

2.633 

50.4006 

1.8532 

76.3632 

2.5664 

9 

Overall  Mean 

499.979 

3.4392 

286.9205 

3.3919 

322.376 

2.6573 

(0,+) 

Mean  STD 

97.8993 

2.2979 

85.2177 

2.4516 

100.2291 

3.0755 

F-11 


rroportio 


APPENDIX  G 

Dispersion  Distribution  Comparisons  Down  Length  of  DZ 


15,000  Ft.  Trail  Distance 
23,500  Ft.  Trail  Distance 


Figure  54.  Dispersion  Distribution  Down  DZ  With  No  Lateral  Separation 
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